-
介绍
-
创建 for 循环函数
-
分配列名
-
方差模型分析
-
设置列名
-
舍入值
-
分配星号
-
在 Excel 工作表中创建表格列表
-
打印最终方差分析对象
-
合并所有变量
-
-
打印方差分析表
介绍
在生物科学中，统计表是数据分析和报告的重要组成部分。用于计算和呈现具有统计意义的发现的传统手动程序既耗时又容易出错。最重要的统计表之一是方差分析 (ANOVA)，可以在广泛的农业、生物和医学研究中找到它。
在这里，我们将重点介绍一种新方法，该方法能够通过几个简单的步骤生成可发布的方差分析表。您可以将此表格导出为 MS word、PDF 和 HTML 文档。如果没有相应编程语言的高级知识，任何统计软件都无法生成此类表格。这将帮助研究人员快速创建方差分析表，同时通过消除手动输入每个命令的需要来节省时间。这比用于计算和表格显示的典型手动过程具有显着优势，后者不仅耗时而且容易出错。
让我们开始吧！
此处列出了在使用某些功能之前必须加载的包。使用library()函数加载这些包。
library(readxl)
library(dplyr)
library(tibble)
library(flextable)数据集中的第一个变量指定重复次数，后面两个是因子变量，其余变量是响应变量。为了评估治疗的相关性，我们将使用双因素随机完全区组设计模型。要导入数据集，我们应该首先使用路径参数定位数据文件。当 col_names 参数设置为 TRUE 时，数据集的第一行将用作变量名。数据文件的元素将作为对象数据存储在 R 中。
data = readxl::read_excel(
path = "Data.xlsx",
col_names = TRUE
)
head(data)# # A tibble: 6 x 10
# Rep Water Priming `Plant height` Spikelets `Spike length` `Grain per spike`
# <dbl> <chr> <chr> <dbl> <dbl> <dbl> <dbl>
# 1 1 CONV NP 100 204 12 35
# 2 1 CONV HP 103 293 15.3 59
# 3 1 CONV OP 95 322 17.5 64
# 4 1 FLOOD NP 99 284 12 34
# 5 1 FLOOD HP 86 304 12.2 44
# 6 1 FLOOD OP 88 326 16.8 69
# # ... with 3 more variables: grain weight <dbl>, Biological yield <dbl>,
# # Grain yield <dbl>当您将 str() 函数应用于对象数据时，您会注意到这两个因子变量在 R 中被识别为字符变量。我们可以使用 as.factor() 函数将这些变量转换为因子。
# Convert categorical variables to factor variables
data$Rep = as.factor(data$Rep)
data$Water = as.factor(data$Water)
data$Priming = as.factor(data$Priming)
str(data)# tibble [27 x 10] (S3: tbl_df/tbl/data.frame)
# $ Rep : Factor w/ 3 levels "1","2","3": 1 1 1 1 1 1 1 1 1 2 ...
# $ Water : Factor w/ 3 levels "CONV","FLOOD",..: 1 1 1 2 2 2 3 3 3 1 ...
# $ Priming : Factor w/ 3 levels "HP","NP","OP": 2 1 3 2 1 3 2 1 3 2 ...
# $ Plant height : num [1:27] 100 103 95 99 86 88 91 89 95 97 ...
# $ Spikelets : num [1:27] 204 293 322 284 304 326 274 291 313 252 ...
# $ Spike length : num [1:27] 12 15.3 17.5 12 12.2 ...
# $ Grain per spike : num [1:27] 35 59 64 34 44 69 32 37 48 32 ...
# $ grain weight : num [1:27] 37 46 56.5 37.5 49.8 59 39.5 44.8 50.6 40.5 ...
# $ Biological yield: num [1:27] 13.3 16.5 17.4 12.3 15 ...
# $ Grain yield : num [1:27] 4.2 5.3 5.85 4.5 5.3 6.1 4.8 5.2 6.3 4.03 ...还将 attach() 函数应用于数据对象以屏蔽其组件。
attach(data)创建 for 循环函数
接下来是困难的部分，因为我们将构建一个 for 循环函数，该函数将询问大量参数，一旦得到回答，将为数据集中的所有响应变量生成一个可发布的 ANOVA 表。
分配列名
首先，我们需要构建一个对象，该对象具有所有数据集响应变量的列名。此数据集中的第一个响应变量位于第四列。因此，响应变量的值将从数据集中的第 4 列分配到第 n 列。
cols <- names(data)[4:ncol(data)]
cols# [1] "Plant height" "Spikelets" "Spike length" "Grain per spike"
# [5] "grain weight" "Biological yield" "Grain yield"方差模型分析
方差分析模型将在 lapply 函数的帮助下实现，该函数返回一个与 X 长度相同的列表，其中包含响应变量的列名。下一步是定义方差分析模型。termlabels 参数中的模型是打印在波浪号后右侧的方差分析模型的组成部分。根据实验设计，您可以提供不同的模型项。
aov.model <- lapply(cols, FUN = function(x)
aov(reformulate(termlabels = "Rep + Water*Priming",
response = x),
data = data))
anova(aov.model[[2]])# Analysis of Variance Table
#
# Response: Spikelets
# Df Sum Sq Mean Sq F value Pr(>F)
# Rep 2 15894.9 7947.5 11.2003 0.0009083 ***
# Water 2 533.1 266.5 0.3756 0.6927466
# Priming 2 28332.7 14166.3 19.9646 4.486e-05 ***
# Water:Priming 4 918.1 229.5 0.3235 0.8581023
# Residuals 16 11353.2 709.6
# ---
# Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1我将选择方差分析表中的第一列、第三列和第五列分别代表自由度、均方和概率值。如果您希望打印 F 值列而不是均方，请使用 4 而不是 3。
# print df, MS and Pvalue
final = anova(aov.model[[2]])[,c(1,3,5)]
final# Df Mean Sq Pr(>F)
# Rep 2 7947.5 0.0009083 ***
# Water 2 266.5 0.6927466
# Priming 2 14166.3 4.486e-05 ***
# Water:Priming 4 229.5 0.8581023
# Residuals 16 709.6
# ---
# Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1我们还需要创建一个对象，该对象应在方差分析表中指定行名，因为我们需要重命名最终发布的方差分析表中的变异源。
# Getting rownames
rnames = rownames(final)
rnames# [1] "Rep" "Water" "Priming" "Water:Priming"
# [5] "Residuals"设置列名
为 ANOVA 表设置与我们已经选择的列相匹配的列名称。
# Setting column names
colnames(final) = c("DF", "MS", "P-value")
colnames(final)[2] = cols[2]
final# DF Spikelets P-value
# Rep 2 7947.5 0.00091
# Water 2 266.5 0.69275
# Priming 2 14166.3 0.00004
# Water:Priming 4 229.5 0.85810
# Residuals 16 709.6舍入值
然后方差分析表中的数据将四舍五入到小数点后三位，输出将保存为数据框。
# Rounding values to 2 decimal place
final = as.data.frame(round(final, digits = 2))
final# DF Spikelets P-value
# Rep 2 7947.47 0.00
# Water 2 266.55 0.69
# Priming 2 14166.35 0.00
# Water:Priming 4 229.53 0.86
# Residuals 16 709.58 NA分配星号
我们不会为每个参数添加 p 值列，而是使用单星号、双星号和 ns 来表示显著性（p<0.05), 高度显著 (p<0.01), 和不重要的 (p>0.05) 影响。之后，我们将合并均方列和我们刚刚分配符号的列。给出列号后，两列将使用粘贴功能合并。
# Assign astericks according to p values
final$sign[final$`P-value` < 0.05] <- "*"
final$sign[final$`P-value` < 0.01] <- "**"
final$sign[final$`P-value` > 0.05] <- "ns"
# Merge MS and significance column together
final[[2]] = paste(final[[2]], ifelse(is.na(final[[4]]), "", final[[4]]))
final# DF Spikelets P-value sign
# Rep 2 7947.47 ** 0.00 **
# Water 2 266.55 ns 0.69 ns
# Priming 2 14166.35 ** 0.00 **
# Water:Priming 4 229.53 ns 0.86 ns
# Residuals 16 709.58 NA <NA>我们将在连接函数之前使用减号从 ANOVA 表中删除额外的列。
final = final[-c(3,4)]
final# DF Spikelets
# Rep 2 7947.47 **
# Water 2 266.55 ns
# Priming 2 14166.35 **
# Water:Priming 4 229.53 ns
# Residuals 16 709.58在 Excel 工作表中创建表格列表
下一个具有挑战性的方面是打印这个最终对象，如果直接打印它会导致单个变量输出而不是所有变量的列表。因此，我决定首先使用 write_excel() 函数为每个变量分别在 excel 表中写入该对象。我通过粘贴每个变量的列名加上带有扩展名的附加文本“-ANOVA”来指定路径参数.xlsx。这将导致在 excel 表中为每个变量打印最终形状的方差分析表。
anova = writexl::write_xlsx(final, path = paste(cols[i], '-ANOVA.xlsx'))打印最终方差分析对象
我能够使用上面创建的包含单独 excel 工作表中的方差分析表列表的对象以发布就绪格式打印最终形状的方差分析表。在打印最终表格之前，我列出了所有带有模式-ANOVA和扩展名的 Excel 表格。.xlsx然后，我使用 lapply 函数通过将 FUN 参数指定为 read_excel 来读取此列表。
# Print final ANOVA object
file.list <- list.files(pattern='*-ANOVA.xlsx')
df.list <- lapply(X = file.list, FUN = read_excel)
df.list# [[1]]
# # A tibble: 5 x 2
# DF `Biological yield`
# <dbl> <chr>
# 1 2 1.93 ns
# 2 2 0.2 ns
# 3 2 46.74 **
# 4 4 0.18 ns
# 5 16 0.7
#
# [[2]]
# # A tibble: 5 x 2
# DF `Grain per spike`
# <dbl> <chr>
# 1 2 24.41 ns
# 2 2 21.12 ns
# 3 2 1374.04 **
# 4 4 30.41 ns
# 5 16 42.05
#
# [[3]]
# # A tibble: 5 x 2
# DF `grain weight`
# <dbl> <chr>
# 1 2 79.42 **
# 2 2 3.37 ns
# 3 2 836.57 **
# 4 4 7 ns
# 5 16 5.23
#
# [[4]]
# # A tibble: 5 x 2
# DF `Grain yield`
# <dbl> <chr>
# 1 2 0.27 **
# 2 2 0.2 *
# 3 2 8.58 **
# 4 4 0.02 ns
# 5 16 0.03
#
# [[5]]
# # A tibble: 5 x 2
# DF `Plant height`
# <dbl> <chr>
# 1 2 0.45 ns
# 2 2 29.93 ns
# 3 2 0.5 ns
# 4 4 15.85 ns
# 5 16 15.5
#
# [[6]]
# # A tibble: 5 x 2
# DF `Spike length`
# <dbl> <chr>
# 1 2 8.71 *
# 2 2 1.1 ns
# 3 2 40.3 **
# 4 4 1.44 ns
# 5 16 1.27
#
# [[7]]
# # A tibble: 5 x 2
# DF Spikelets
# <dbl> <chr>
# 1 2 7947.47 **
# 2 2 266.55 ns
# 3 2 14166.35 **
# 4 4 229.53 ns
# 5 16 709.58合并所有变量
在接下来的阶段，我使用了 list.cbind() 函数来连接之前生成的所有列表。此组合方差分析表包含数据集中所有响应变量的 DF 和 MS 列。但是，我们需要删除每个变量的 DF 列。在以下命令中，我创建了一个包含重复列名称的对象，随后通过在该对象前使用减号删除了重复列。
# Combined ANOVA table for all variables
aov.table = rlist::list.cbind(df.list)
# Remove duplicate columns for DF
dup.cols = which(duplicated(names(aov.table)))
aov.table = aov.table[,-dup.cols]
aov.table# DF Biological yield Grain per spike grain weight Grain yield Plant height
# 1 2 1.93 ns 24.41 ns 79.42 ** 0.27 ** 0.45 ns
# 2 2 0.2 ns 21.12 ns 3.37 ns 0.2 * 29.93 ns
# 3 2 46.74 ** 1374.04 ** 836.57 ** 8.58 ** 0.5 ns
# 4 4 0.18 ns 30.41 ns 7 ns 0.02 ns 15.85 ns
# 5 16 0.7 42.05 5.23 0.03 15.5
# Spike length Spikelets
# 1 8.71 * 7947.47 **
# 2 1.1 ns 266.55 ns
# 3 40.3 ** 14166.35 **
# 4 1.44 ns 229.53 ns
# 5 1.27 709.58在最后阶段，我们只需要在方差分析表中添加指定变异来源的行名。
# Names for sources of variation in ANOVA
rownames(aov.table) = rnames
aov.table# DF Biological yield Grain per spike grain weight Grain yield
# Rep 2 1.93 ns 24.41 ns 79.42 ** 0.27 **
# Water 2 0.2 ns 21.12 ns 3.37 ns 0.2 *
# Priming 2 46.74 ** 1374.04 ** 836.57 ** 8.58 **
# Water:Priming 4 0.18 ns 30.41 ns 7 ns 0.02 ns
# Residuals 16 0.7 42.05 5.23 0.03
# Plant height Spike length Spikelets
# Rep 0.45 ns 8.71 * 7947.47 **
# Water 29.93 ns 1.1 ns 266.55 ns
# Priming 0.5 ns 40.3 ** 14166.35 **
# Water:Priming 15.85 ns 1.44 ns 229.53 ns
# Residuals 15.5 1.27 709.58当我完成我的代码时，我需要编写一个for 循环函数来为数据集中的所有响应变量重复它。这可以通过将前面的所有指令包含在 for 循环函数中来实现，我在其中为 i 在数据集的 1 到 n 列中指定了一个范围，前三个变量除外。
打印方差分析表
最后，在 MS Word 文档中，我们必须打印此可发布的方差分析表。为此，我们将利用 flextable 包中的 flextable 函数。我们将在数据参数中提供 aov.table 对象，我们还可以在管道运算符之后为表的行名提供列名。通过在 bold() 函数中使用 part 参数作为标题，您可以将列名加粗。
单击编译报告图标以多种文档格式生成此表。您可以从多种文档类型中进行选择。我要选择 MS Word 文档。
这就是方差分析表最终的样子。您认为经过如此繁琐的工作后，只需点击几下鼠标即可在 R 中生成可供发布的方差分析表吗？答案是肯定的，您只需为自己的数据更改几个参数即可获得相同的输出。
当我们将变量从字符转换为因子时，您必须在开头部分提供数据集中因子变量的名称。然后您需要选择一个最能描述实验设计的模型项。以下是不同实验设计的一些模型项。您可以选择最能描述您的实验设计的一个。您必须更改的另一件事是选择可能在您的数据集中不同的列名称的适当范围。就这样！
扫一扫
331
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
