7 第七讲:如何进行基本的数据分析: t-test和ANOVA

{r setup, include=FALSE}2 knitr::opts_chunk$set(echo = TRUE) library(tidyverse) # Wickham的数据整理的整套工具 pdf.options(height=10/2.54, width=10/2.54, family="GB1") # 注意:此设置要放在最后 ## 语法实现

7.0.1 t检验

7.0.1.1 理论基础

我们来讲讲t检验。最简单的是单样本t检验,它将样本的均值与某个固定值进行比较,这个固定的值通常被认为是该样本所在的总体。在这种情境下,按照NHST的逻辑,我们首先假设\(H_0\)为真,即样本确实是从总体均值为这个固定值的分布中抽样出来的。然后根据\(H_0\)的统计模型来计算当前数据或者更加极端数据的概率。如果概率很低,那么在一次抽样中这样的小概率事件不太可能发生,因此拒绝假设\(H_0\)\(p\)值与假设水平进行比较时,我们假设\(H_0\)为真。我们有三个选择,等于、大于或小于。如果不等于,那么是双样本检验,如果小于或大于,则为单样本检验。

下面具体来介绍这些检验的不同类型。 首先是单样本检验,当p值小于等于\(\mu\)时,\(H_0\)大于等于\(\mu\)。另一种单样本t检验是当\(H_0\)大于等于\(\mu\)时。

接下来是双样本检验,其零假设是:两个样本来自同一个总体。我们可以用p值或t值检验来比较两个样本之间的差异。如果两个样本的数据是高度相关的,使用配对样本t检验。如果两个样本的数据是毫不相关即完全独立的,则使用独立样本t检验。

在R语言中,有多个函数可以用于t检验,这些函数包括在R-base包,也包括在其他的包中。我们将会介绍Bruce R包,它是为心理学数据分析而进行优化过的,更方便。

Bruce R包集成了很多底层R包,有机整合了其中的R函数,功能主要覆盖: 1.基础R编程(自动设置文件夹路径、一站式数据导入导出、数据匹配拼接等) 2.多变量计算(极简化计算多题项平均分/总分及反向计分、数值重新编码等) 3.信度与因素分析(量表信度分析、探索性因素分析EFA、主成分分析PCA、验证性因素分析CFA) 4.描述与相关分析(描述统计、频数统计、相关分析、相关系数差异性检验) 5.t检验与方差分析(单样本/独立样本/配对样本t检验、多因素被试间/被试内/混合设计方差分析ANOVA、简单效应检验与多重比较) 6.普通与多层回归分析(多种回归模型结果的完整报告与表格输出、多层线性模型HLM补充分析) 7.中介与调节效应分析(普通与多水平的中介效应、调节效应、有调节的中介效应等,包括简单斜率分析、条件中介作用、链式中介作用、跨层调节作用等) 8.统计与绘图辅助工具(总均值/组均值中心化处理、时间序列交叉相关分析与格兰杰因果检验、ggplot2绘图主题theme_bruce、颜色卡/配色方案展示等)

####解决问题 问题: 1. “在penguin数据中,参与者对于家的依恋水平是否小于/等于/大于均值水平(假设在总体水平上,人们对家庭的依恋水平(HOME)均值为3.5)?” 2. “在penguin数据中,男女生在亲密关系经验(ECR)中的得分是否存在显著差异?” 3. “在penguin数据中,参与者在不同时间的温度差异是否具有显著性?”

library(bruceR)
## 
## bruceR (v2023.9)
## Broadly Useful Convenient and Efficient R functions
## 
## Packages also loaded:
## ✔ data.table ✔ emmeans
## ✔ dplyr      ✔ lmerTest
## ✔ tidyr      ✔ effectsize
## ✔ stringr    ✔ performance
## ✔ ggplot2    ✔ interactions
## 
## Main functions of `bruceR`:
## cc()             Describe()  TTEST()
## add()            Freq()      MANOVA()
## .mean()          Corr()      EMMEANS()
## set.wd()         Alpha()     PROCESS()
## import()         EFA()       model_summary()
## print_table()    CFA()       lavaan_summary()
## 
## For full functionality, please install all dependencies:
## install.packages("bruceR", dep=TRUE)
## 
## Online documentation:
## https://psychbruce.github.io/bruceR
## 
## To use this package in publications, please cite:
## Bao, H.-W.-S. (2023). bruceR: Broadly useful convenient and efficient R functions (Version 2023.9) [Computer software]. https://CRAN.R-project.org/package=bruceR
## 
## These packages are dependencies of `bruceR` but not installed:
## - ggtext, vars, phia, GPArotation
## 
## ***** Install all dependencies *****
## install.packages("bruceR", dep=TRUE)
library(dplyr)
library(tidyr)
df.pg.raw <-  read.csv('./data/penguin/penguin_rawdata.csv',
                       header = T, sep=",", stringsAsFactors = FALSE)

首先载入所需的包。BruceR包载入后会有一大串提示信息,包括载入的包、主要函数和网址等。与先前操作相似,我们使用相对路径读取数据,可以查看前5行。

df.pg.mean <- df.pg.raw %>%
  dplyr::filter(sex > 0 & sex < 3) %>%
  dplyr::mutate(ECR_mean = rowMeans(select(., starts_with("ECR")),na.rm = T),
                HOME_mean = rowMeans(select(., starts_with("HOME")),na.rm = T),
                sex=as.factor(sex)
                ) %>%
  dplyr::select(sex, ECR_mean, HOME_mean, Temperature_t1,Temperature_t2)

其次,我们需要进行数据预处理,计算家庭依恋的均值。这需要使用之前学过的数据预处理函数和方法,我们将使用管道进行操作。为了计算第二个问题,我们首先筛选出两个性别,即男性和女性。然后,我们对数据进行操作,计算两个变量的得分均值,用mutate生成两个新变量。我们使用命令na.rm来去除包含缺失值的数据,以避免出现错误。我们可以使用命令将性别转换为factor,以便在进行t-test时更容易分类。在使用R语言时,不同的函数可能会有微小的区别,因此我们需要注意这些细节。在选择数据时,应该只选取感兴趣的变量进行分析,如年龄、性别、亲密关系、问卷得分和温度测量值等。

进行t检验时,需要输入数据框、y变量和x变量,同时需要注意是否是配对样本和是否满足方向性。

问题一(单样本t检验):

  bruceR::TTEST(df.pg.mean, "HOME_mean", test.value = 3.5, test.sided = ">")
## 
## One-Sample t-test
## 
## Hypothesis: one-sided (μ > 3.5)
## t approximation invoked.
## Descriptives:
## ───────────────────────────
##   Variable    N Mean (S.D.)
## ───────────────────────────
##  HOME_mean 1490 3.99 (0.73)
## ───────────────────────────
## 
## Results of t-test:
## ───────────────────────────────────────────────────────────────────────────────────────────────────
##                                   t   df     p     Difference [95% CI] Cohen’s d [95% CI]      BF10
## ───────────────────────────────────────────────────────────────────────────────────────────────────
## HOME_mean: (HOME_mean - 3.5)  25.97 1489 <.001 ***    0.49 [0.46, Inf]   0.67 [0.63, Inf] 2.42e+119
## ───────────────────────────────────────────────────────────────────────────────────────────────────
  bruceR::TTEST(df.pg.mean, "HOME_mean", test.value = 3.5, test.sided = "<")
## 
## One-Sample t-test
## 
## Hypothesis: one-sided (μ < 3.5)
## t approximation invoked.
## Descriptives:
## ───────────────────────────
##   Variable    N Mean (S.D.)
## ───────────────────────────
##  HOME_mean 1490 3.99 (0.73)
## ───────────────────────────
## 
## Results of t-test:
## ──────────────────────────────────────────────────────────────────────────────────────────────────
##                                   t   df     p     Difference [95% CI] Cohen’s d [95% CI]     BF10
## ──────────────────────────────────────────────────────────────────────────────────────────────────
## HOME_mean: (HOME_mean - 3.5)  25.97 1489 1.000       0.49 [-Inf, 0.52]  0.67 [-Inf, 0.72] 6.84e-04
## ──────────────────────────────────────────────────────────────────────────────────────────────────
  result.ttest <- capture.output({
  bruceR::TTEST(data=df.pg.mean, 
                y="HOME_mean", 
                test.value = 3.5, 
                test.sided = "=",
                file = "./output/chp7/single_t.doc")
})
writeLines(result.ttest, "./output/chp7/single_t.md") # .md最整齐

mean difference代表两组数据的均值差异,Cohen's *d*是t-test的一个效应量。test side表示双尾还是单尾。test value是要检验的变量。factor variance表示是否对因素进行反转,digital参数指小数点后保留的位数,file参数用于保存结果到word文档。

使用capture output函数可以整理结果并复制到results.ttest变量中。将结果复制到名为results.ttest的变量中,并将结果写入一个名为single_t.md的文件中。在根目录中打开R4Psy项目文件,我们可以看到已经运行过的输出结果,其中包括一个符合心理学格式的三线表格,其中包含了均值、t值、自由度、p值、协方差和效应量等信息。这个表格非常干净。

问题二(独立样本t检验):

    result.ttest <- capture.output({   
    bruceR::TTEST(data=df.pg.mean, 
                y="ECR_mean",
                x="sex")
})

writeLines(result.ttest, "./output/chp7/inde_t.md")

然后是独立样本t检验,比较男性和女性之间的差异。使用了ECR_mean作为样本得分,用性别作为分组变量。结果是two-side的,有描述性的结果,结果包括每组的n、均值和方差齐性检验。如果方差不齐性,可以使用校正方法。t检验的结果,包括t值、df值、p值、mean difference、Cohen’s d值。

问题三(配对样本t检验):

   result.ttest <- capture.output({
   bruceR::TTEST(data=df.pg.mean, 
                 y = c("Temperature_t1",              
                       "Temperature_t2"),
                 paired = T)  #是否为配对样本t检验?默认是FALSE
})
writeLines(result.ttest, "./output/chp7/pair_t.md")

配对样本t检验也很简单,只需要将y变成两个值,然后使用t1和t2进行检验。

请大家注意,我们之前讲过如何管理和放置你的R-Project,可以避免路径错误。如果你不熟悉路径,可以从github上下载整个文件夹,打开R4psy文件夹中的R-Project文件。我们的助教已经在群里分享了GitHub链接。打开链接后,点击绿色按钮下面的“code”,再点击“download zip”,就可以下载整个文件夹。下载后解压,把文件夹放在一个完整的文件夹里。文件夹里有一个叫做“R4Psy.Rproj”的文件,双击它就可以打开Rstudio。

在console里输入代码时,如果按回车后,有时可能看到一个大于号“>”或加号”+“。其中“+”表示代码没有输完,可能是缺少反括号或者某些代码没有输入完整。可以补上反括号或按退出键(“ESC”,一般在键盘的左上角)退出。Rstudio是一个好的IDE,可以帮助我们导航代码。在Rmarkdown里,代码块以chunk1、chunk2等命名,可以快速定位到代码框。有些chunk没有命名,只有编号。也可以通过outline来定位。

如果你想要快速运行代码,推荐使用快捷键control+回车,在Mac上是command+回车。也可以通过查找keyboard shortcut来找到更多快捷键。

在进行统计分析后,可以使用软件包reports来标准化报告结果,包括均值、标准差、t值、p值和Cohen’s d等。

7.0.2 方差分析

我们将使用BruceRTidyverse软件包,但需要注意,安装BayesFactor包可能会出现错误,因为它依赖于Rtools包。如果出现错误,可以在官方网站上搜索如何安装BayesFactor包。注意安装依赖包,否则可能会出现问题。在使用BruceR时,可能会遇到一些警告和问题,需要自己解决。

7.0.2.1 拟解决问题

  1. “在match数据中,被试在匹配任务和不匹配任务上的反应时是否存在差异?”
  2. “在match数据中,被试在不同身份(self vs other)与不同效价(moral vs immoral)的条件组合下反应时是否存在差异?”
df.match <-  read.csv('./data/match/match_raw.csv',
                      header = T, sep=",", stringsAsFactors = FALSE) %>%
  # 拆分单元格内字符串
  tidyr::separate(col=Shape,
                  into=c("Valence","Identity"),
                  sep="(?<=moral|immoral)(?=Self|Other)") %>% 
  dplyr::select(Sub, Valence, Identity, everything()) %>%
  dplyr::filter(ACC == 1) %>%    # 只选择回答正确的数据
  dplyr::filter(!is.na(RT)) %>%  # 剔除缺失值
  # remove outliers below and above 3rd sd
  dplyr::filter(RT > quantile(RT, 0.0015) & RT < quantile(RT, 0.9985)) %>%
  dplyr::mutate(RT = as.numeric(RT)) %>%
  dplyr::mutate(Valence = as.factor(Valence),
                Identity = as.factor(Identity))
df.match.raw <-  read.csv('./data/match/match_raw.csv',
                       header = T, sep=",", stringsAsFactors = FALSE)%>%
                 tidyr::separate(., col=Shape,into = c("Valence","Identity"),
                                              sep = "(?<=moral|immoral)(?=Self|Other)")%>% # 拆分单元格内字符串
                 dplyr::select(Sub,Valence,Identity,everything())%>%
                 dplyr::filter(ACC == 1) %>%  # 只选择回答正确的数据
                 dplyr::filter(RT > quantile(RT, 0.0015) & RT < quantile(RT, 0.9985)) %>% # remove outliers below and above 3rd sd
                 dplyr::mutate(RT = as.numeric(RT)) %>%
                 dplyr::mutate(across("ACC", as.factor)) %>%
                 dplyr::mutate(Valence = as.factor(Valence)) %>%
                 dplyr::mutate(Identity = as.factor(Identity)) %>%
                 dplyr::filter(!is.na(RT)) #剔除缺失值
DT::datatable(head(df.match.raw, 10),
              fillContainer = TRUE, options = list(pageLength = 5))

在数据清理时,可以先读取数据,然后进行管道操作,将Shape列分割成Valence和Identity两列,并选择Base、Valence、Identity和Everything列。然后进行数据清理,只选择反应正确的反应时间,并去掉缺失值和极端反应时间。我们使用了三个标准差的方法来去除极端值。我们将反应时间转换为数值型,将情感和身份转换为因子型。我们将所有操作应用于数据框df.match,这是经过预处理的数据。

df.match.mean <- df.match.raw %>%
    dplyr::group_by(Sub, Match) %>%
    dplyr::summarise(
        n = n(),
        rt_sd = sd(as.numeric(RT), na.rm = T),
        rt_mean = mean(as.numeric(RT), na.rm = T),
        acc_mean=mean(as.numeric(as.character(ACC)),na.rm=T)#对于factor数据,先转成character,再变成numeric
        ) %>%
    dplyr::ungroup()
## `summarise()` has grouped output by 'Sub'. You can override using the `.groups`
## argument.

我们需要求出每个被试在匹配和不匹配任务上的均值,才能进行比较。为此,我们使用分组进行预处理,分组的目标是为了得到每个组的描述性统计,如RT的mean和SD,以及正确试次的个数。我们可以使用summarize快速得出这些统计信息。接着,我们可以使用t-test或方差分析来比较匹配和不匹配反应是否有差异。

问题一(单因素被试内设计)之单因素方差分析:

#比较被试在匹配任务和不匹配任务上的反应时是否存在差异
#单因素被试内设计(长型数据)
result.anova <- capture.output({
df_within <- df.match.mean%>%
  bruceR::MANOVA(subID="Sub",#被试id
                 dv="rt_mean",#因变量,因为是MANOVA,实际上因变量可以再后面设置很多个
                 within="Match") #设置条件,因素分析的条件
})
## 
##     * Data are aggregated to mean (across items/trials)
##     if there are >=2 observations per subject and cell.
##     You may use Linear Mixed Model to analyze the data,
##     e.g., with subjects and items as level-2 clusters.
writeLines(result.anova, "./output/chp7/anova_1.md")
# 若不符合球形假设要加上:sph.correction = "GG"

在方差分析中,我们可以使用BruceR包中的MANOVA()函数。参数分别是subID,dv和within,即可直接对其进行操作并输出结果。这几个参数符合心理学的命名规范。我们可以得到match和mismatch两种条件的描述性统计,包括f值、p值、eta2、pash2、95%置信区间和generalized eta2。因为这是一个完全被试内设计的分析,所以不需要进行方差齐性检验。

问题二(两因素被试内设计)之重复测量方差分析:

df.match.mean <- df.match %>%
    dplyr::group_by(Sub,Sex, Valence,Identity) %>%
    dplyr::summarise(
        n = n(),
        rt_sd = sd(as.numeric(RT), na.rm = T),
        rt_mean = mean(as.numeric(RT), na.rm = T)
        ) %>%
    dplyr::ungroup()
## `summarise()` has grouped output by 'Sub', 'Sex', 'Valence'. You can override
## using the `.groups` argument.

对于被试身份和效价这两个因素,我们可以直接从match中选出数值对,然后进行groupby和summarize操作。需要注意的是,我们现在是以subject、variance和identity为基础进行分析。

df.match.within <- df.match.mean %>%
  dplyr::select(-c(n, rt_sd)) %>%
  dplyr::mutate(Valence = paste("A_", Valence, sep = ""),#将变量的名称进行修改转换,不生产新的变量
                Identity = paste("B_", Identity, sep = "")) %>%
  # 将morality和identity组合名称起来生成一个新的变量"Conds"
  tidyr::unite("Conds", Valence:Identity, sep = "&",remove=TRUE) %>% 
  tidyr::pivot_wider(names_from = Conds,
                     values_from = rt_mean)

head(df.match.within)
## # A tibble: 6 × 6
##     Sub Sex    `A_immoral&B_Other` `A_immoral&B_Self` `A_moral&B_Other`
##   <int> <chr>                <dbl>              <dbl>             <dbl>
## 1  7302 female               0.706              0.724             0.662
## 2  7303 male                 0.741              0.754             0.720
## 3  7304 female               0.756              0.730             0.757
## 4  7305 male                 0.685              0.659             0.641
## 5  7306 male                 0.778              0.763             0.660
## 6  7307 female               0.751              0.677             0.716
## # ℹ 1 more variable: `A_moral&B_Self` <dbl>
#str(df.match.within)
#被试在不同身份(self vs other)与不同效价(moral vs moral)的条件组合下反应时是否存在差异)
  result.anova <- capture.output({
  res_rmANOVA_1 <- bruceR::MANOVA(data=df.match.within,
                   dvs="A_immoral&B_Other:A_moral&B_Self",
                   dvs.pattern="A_(.+)B_(.+)",#正则表达式
                   within=c("A_","B_"))#表示2个主条件
  })
## 
## Note:
## dvs="A_immoral&B_Other:A_moral&B_Self" is matched to variables:
## A_immoral&B_Other, A_immoral&B_Self, A_moral&B_Other, A_moral&B_Self
writeLines(result.anova, "./output/chp7/anova_2.md")

如果我们有的是宽数据,我们可以先将宽数据转成长数据,然后选择变量并按照Bruce R的命名方式在前面加上a和b表示自变量a和自变量b。接着,我们可以使用unite函数将Valence和Identity合并成一个新的变量,再使用pivot wide将其转换为更宽的形式。这样,我们就可以看到每个条件下的性别、被试id和条件。 最后,我们将结果写成md格式,得到一个简洁的表格,其中包含F值、p值、partial eta square、generalized square和partial eta square的95%置信区间。长数据也可以使用同样的方法处理。

#被试在不同身份(self vs other)与不同效价(good vs bad)的条件组合下反应时是否存在差异)
#head(df.match.mean)
  result.anova <- capture.output({
  res_rmANOVA_2 <- bruceR::MANOVA(data=df.match.mean,
                   dv="rt_mean",
                   within=c("Valence","Identity"),
                   subID="Sub")
  })
## 
##     * Data are aggregated to mean (across items/trials)
##     if there are >=2 observations per subject and cell.
##     You may use Linear Mixed Model to analyze the data,
##     e.g., with subjects and items as level-2 clusters.
writeLines(result.anova, "./output/chp7/anova_3.md")

简单效应分析:

result.check <- capture.output({
sim_eff_1 <- res_rmANOVA_1 %>%
  bruceR::EMMEANS("A_", by="B_") # 简单效应分析
})
writeLines(result.check, "./output/chp7/check.md")

#sim_eff_1 <- res_rmANOVA_1 %>%
  #EMMEANS("B_", by="A_")   和上一个同理,只是转换了不同的形式

除了基础分析方法,我们还可以使用EMMEANS进行简单效应分析和多重比较矫正等操作。通过输入之前的方差分析结果,我们可以使用EMMEANS a by b来进行简单效应分析,以比较在不同条件下的主效应是否存在差异。

我们要分析a在b的不同条件下的效应,包括other和self条件。在other条件下,a的效应是b条件的一个主效应,我们称之为简单效应。我们可以看到,在b的两个条件下,a的简单效应实际上应该有一个交互作用。在other条件下,它不显著,在self条件下,它很显著。 此外,我们还可以看到a的quantity,即在self条件下,a的主效应moral和immoral之间差异非常大,quantity大于0.8,是一个大的效应。