用于计算平均值和标准误差的笨拙代码答案

【问题标题】：Clumpsy code for calculating mean and standard error用于计算平均值和标准误差的笨拙代码
【发布时间】：2020-12-08 15:21:28
【问题描述】：

我试图编写一段代码来计算我的数据的平均值和标准误差并将其放入一个新的小标题中。

但是，感觉非常笨拙。有谁知道可以使我的代码更优雅的包或其他技巧？

我需要计算多个子组 (days_incubated) 的均值和 se。

library(dplyr)
library(tibble)
library(tidyr)
library(data.table)

library(plotrix)

df2 <- df1%>%
  group_by(days_incubated)%>%
  summarise_each(funs(mean, se= std.error))%>% # Calculating mean and standard error
  mutate_if(is.numeric, round, digits = 2) # Round off the data

df2_trans <- transpose(df2) # Transposing data table

colnames(df2_trans) <- rownames(df2) # Get row and colnames in order
rownames(df2_trans) <- colnames(df2) # Get row and colnames in order                

df2_trans <- rownames_to_column(df2_trans, "mass") # Making row names into a column 

df3_trans <- df2_trans%>% # Converting one column into two
  separate(mass, c("mass","type"), sep = "([_])")

mean_target <- c("mean", "incubated")

mean <- df3_trans%>% # Mean table
  filter(type %in% mean_target)%>%
  rename("mean day 0"="1")%>%
  rename("mean day 4"="2")%>%
  rename("mean day 10"="3")%>%
  rename("mean day 17"="4")%>%
  rename("mean day 24"="5")%>%
  rename("mean day 66"="6")%>%
  rename("mean day 81"="7")%>%
  rename("mean day 94"="8")%>%
  rename("mean day 116"="9")%>%
  select("mass", "mean day 0", "mean day 4", "mean day 10", "mean day 17", "mean day 24", "mean day 66", "mean day 81", "mean day 94", "mean day 116")%>%
  slice(-c(1))

se_target <- c("se", "incubated")

se <- df3_trans%>% # SE table
  filter(type %in% se_target)%>%
  rename("se day 0"="1")%>%
  rename("se day 4"="2")%>%
  rename("se day 10"="3")%>%
  rename("se day 17"="4")%>%
  rename("se day 24"="5")%>%
  rename("se day 66"="6")%>%
  rename("se day 81"="7")%>%
  rename("se day 94"="8")%>%
  rename("se day 116"="9")%>%
  select("mass", "se day 0", "se day 4", "se day 10", "se day 17", "se day 24", "se day 66", "se day 81", "se day 94", "se day 116")%>%
  slice(-c(1))

# join mean and se tables
mean_se <- mean %>% #merging mean and se dataset
  full_join(se, by=("mass"))%>%
  select("mass","mean day 0","se day 0", "mean day 4", "se day 4", "mean day 10", "se day 10", "mean day 17", "se day 17", "mean day 24", "se day 24", "mean day 66", "se day 66", "mean day 81", "se day 81", "mean day 94", "se day 94", "mean day 116", "se day 116") # Putting columns in correct order

这是数据：

df1 <- structure(list(days_incubated = c("0", "0", "0", "0", "0", "4", 
"4", "4", "4", "4", "10", "10", "10", "10", "10", "17", "17", 
"17", "17", "17", "24", "24", "24", "24", "24", "66", "66", "66", 
"66", "66", "81", "81", "81", "81", "81", "94", "94", "94", "94", 
"94", "116", "116", "116", "116", "116"), i.x33.031 = c(7.45, 
0, 78.2, 16.49, 18.77, 104.5, 28.95, 26.05, 4.11, 62.09, 1.95, 
6.75, 1.41, 3.34, 3.02, 0, 100.28, 0.2, 32.66, 0, 0, 370.57, 
7.24, 133.63, 55.26, 0.16, 5.5, 25.17, 16.59, 3.3, 23.95, 30.61, 
4.04, 0, 6.58, 0.08, 0.01, 0, 0.38, 0, 0, 0, 0, 0.18, 0), i.x35.034 = c(0, 
0, 0.15, 0.02, 0.01, 0.04, 0.04, 0.05, 0.02, 0.09, 0.02, 0, 0.04, 
0.01, 0, 0, 0.22, 0, 0.08, 0, 0, 0.66, 0.01, 0.2, 0.12, 0.01, 
0.01, 0.04, 0.01, 0.01, 0.01, 0.04, 0, 0, 0, 0, 0, 0, 0.01, 0, 
0, 0.02, 0, 0, 0.02), i.x36.017 = c(0.47, 0.09, 0.28, 0.02, 0.03, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.05, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.3, 0.06, 0.32, 0, 0, 0, 0, 0.12, 
0, 0.02), i.x39.959 = c(0.02, 0, 0.08, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.04, 0, 0, 0, 0, 0, 0.01, 0, 0, 
0, 0, 0, 0.01, 0.02, 0.06, 0.03, 0.03, 0, 0, 0.02, 0.01, 0, 0, 
0), i.x40.023 = c(0.35, 0.02, 0.48, 0.06, 0, 1.25, 0.09, 0.1, 
0.03, 0, 0.09, 0.07, 0.55, 0.09, 0.07, 0, 0.63, 0, 0.09, 0.07, 
0.02, 1.11, 0.04, 0.59, 0.13, 0, 0.01, 0.02, 0, 0, 0, 0, 0.01, 
0.02, 0.06, 0.01, 0.01, 0.01, 0.01, 0.04, 0, 0.08, 0, 0, 0.01
)), row.names = c(NA, -45L), class = "data.frame")

【问题讨论】：

这段代码对我来说不起作用。您是否使用来自data.table 的transpose()？否则我假设它来自purrr，因为你加载了tidyverse。
我编辑了这些包，并且代码似乎按照您的预期运行。否则，请随时将其编辑回来。

标签： r

【解决方案1】：

这能满足您的需求吗？它似乎重现了输出。请注意，我认为您的结果实际上标有错误的列。您的其中一个步骤对列进行了字符排序，因此它们与您稍后重命名时所期望的顺序不同。

我们首先使用pivot_longer()，这样我们就有了整洁的数据。从那里，我们可以分组和总结以计算平均误差和标准误差。然后我们pivot_wider() 将其移回您的结果的宽格式。

library(dplyr)
library(tidyr)
library(plotrix)

col_order <- paste0(rep(paste0(c("mean", "se"), "_day_"), length(unique(df1$days_incubated))),
                    rep(unique(df1$days_incubated), each = 2))

df1 %>% 
  pivot_longer(-days_incubated, names_to = "mass") %>% 
  group_by(days_incubated, mass) %>% 
  summarize(mean = mean(value),
            se = std.error(value), .groups = "drop") %>%
  pivot_wider(names_from = days_incubated, values_from = c("mean", "se"),
              names_glue = "{.value}_day_{days_incubated}") %>% 
  relocate(mass, all_of(col_order))

# A tibble: 5 x 19
  mass  mean_day_0 se_day_0 mean_day_4 se_day_4 mean_day_10 se_day_10 mean_day_17 se_day_17 mean_day_24 se_day_24 mean_day_66 se_day_66 mean_day_81
  <chr>      <dbl>    <dbl>      <dbl>    <dbl>       <dbl>     <dbl>       <dbl>     <dbl>       <dbl>     <dbl>       <dbl>     <dbl>       <dbl>
1 i.x3~     24.2    13.9        45.1    17.5          3.29    0.932        26.6     19.5        113.       68.6        10.1       4.67       13.0  
2 i.x3~      0.036   0.0287      0.048   0.0116       0.014   0.00748       0.06     0.0429       0.198     0.121       0.016     0.006       0.01 
3 i.x3~      0.178   0.0867      0       0            0       0             0        0            0.01      0.01        0         0           0    
4 i.x3~      0.02    0.0155      0       0            0       0             0        0            0.008     0.008       0.002     0.002       0.006
5 i.x4~      0.182   0.0978      0.294   0.240        0.174   0.0941        0.158    0.119        0.378     0.210       0.006     0.004       0.018
# ... with 5 more variables: se_day_81 <dbl>, mean_day_94 <dbl>, se_day_94 <dbl>, mean_day_116 <dbl>, se_day_116 <dbl>

逐行

首先，我添加了一个简短的步骤来帮助列排序。

col_order 是我们根据我们正在寻找的正确最终顺序创建的字符向量。基本上我们只是连接字符串以匹配列名。您可以跳过此步骤，只需像在 select() 语句中一样单独输入每个名称，但这样可以节省输入时间。

现在进入管道。

pivot_longer() 将从列中获取数据并将它们堆叠起来。列名现在将位于名为mass 的变量中。默认情况下，这些值将位于名为 value 的新变量中。
group_by() 对指定变量进行分组，这些变量是您希望在最终表格中出现的行和列组合。
summarize() 将每个 mass 和 days_incubated 组合的多行折叠成一个新行，其中包含两列：均值和 se。 .groups = "drop" 参数表示要删除组（即，tibble 不再分组，请查看帮助文件了解更多详细信息）。
pivot_wider() 采用长格式摘要标题并将数据转出到列中，这与您的格式非常相似。它说从days_incubated 中获取列名，并从mean 和se 列中获取值。 names_glue 参数指定如何命名新列。它使用glue 包中的语法，{ 括号基本上表示“在此处替换值”。因此，每个名称都是 '[name of value column]_day_[days_incubated]'。
relocate() 只是对列重新排序。这对于dplyr 来说相对较新，并且是现在重新排列列而不是select() 的首选方式。它表示首先获取mass，然后是all_of() 我们之前创建的col_order 向量中的列。

【讨论】：

感谢@Adam，它看起来很棒！列的顺序非常重要。我将如何使用relocate()
我想知道groups = "drop") 和names_glue = "{.value}_day_{days_incubated}") 是做什么的？
哇，这真的很优雅！
如果您还没有，请查看一些关于“整洁数据”的文献 (cran.r-project.org/web/packages/tidyr/vignettes/tidy-data.html)。当您能够识别数据的结构以及如何改进它时，计算通常会变得更加精简。
没有理由对我的代码做任何事情。在您的原始代码中，查看df2_trans。看第二排。看到列的顺序了吗？它是 0、10、116、17 等等。当您稍后使用rename() 时，您说第一列是第 0 天，第二列是第 4 天，依此类推。这是不正确的，因为数据是在某个地方排序的。

【解决方案2】：

我想提供一个主要在基础 R 中的解决方案。

library(dplyr)

# Define a function that computes the standard error of a vector.
f.SE <- function(x) sd(x)/sqrt(length(x))

# Define a function that compute the standard error of a column (MARGIN = 2 for column). 
f.colSE <- function(mydf) apply(mydf, MARGIN = 2, FUN = f.SE)

# Define a function that combines the means and the standard errors of a column
f.colMeanSE <- function(mycol) cbind(means = colMeans(mycol), se = f.colSE(mycol))

# Apply the function to all numeric columns in your data frame by using the days_incubated subgroup as the indices. The result is a list. 
my_mean_se1 <- by(data = df1[,-1], INDICES = df1$days_incubated, FUN= f.colMeanSE); my_mean_se1

#df1$days_incubated: 0
#means          se
#i.x33.031 24.182 13.91245823
#i.x35.034  0.036  0.02874022
#i.x36.017  0.178  0.08668333
#i.x39.959  0.020  0.01549193
#i.x40.023  0.182  0.09779571
#------------------------------------------------------------ 
#  df1$days_incubated: 10
#means          se
#i.x33.031 3.294 0.932119091
#i.x35.034 0.014 0.007483315
#i.x36.017 0.000 0.000000000
#i.x39.959 0.000 0.000000000
#i.x40.023 0.174 0.094106323
#------------------------------------------------------------

# Bind the columns. 
my_mean_se2 <- do.call(cbind, my_mean_se1); my_mean_se2

#           means          se means          se means          se  means          se   means         se  means          se  means       se  means
#i.x33.031 24.182 13.91245823 3.294 0.932119091 0.036 0.036000000 26.628 19.46476571 113.340 68.5767486 45.140 17.49416646 10.144 4.665935 13.036
#i.x35.034  0.036  0.02874022 0.014 0.007483315 0.008 0.004898979  0.060  0.04289522   0.198  0.1212601  0.048  0.01157584  0.016 0.006000  0.010
#i.x36.017  0.178  0.08668333 0.000 0.000000000 0.028 0.023323808  0.000  0.00000000   0.010  0.0100000  0.000  0.00000000  0.000 0.000000  0.000
#i.x39.959  0.020  0.01549193 0.000 0.000000000 0.006 0.004000000  0.000  0.00000000   0.008  0.0080000  0.000  0.00000000  0.002 0.002000  0.006
#i.x40.023  0.182  0.09779571 0.174 0.094106323 0.018 0.015620499  0.158  0.11939012   0.378  0.2103188  0.294  0.23972067  0.006 0.004000  0.018
 #                  se means         se
#i.x33.031 6.002108463 0.094 0.07304793
#i.x35.034 0.007745967 0.002 0.00200000
#i.x36.017 0.000000000 0.136 0.07194442
#i.x39.959 0.004000000 0.024 0.01122497
#i.x40.023 0.011135529 0.016 0.00600000

# Bind the names of the columns to the names of days_incubated
names_days <- unique(df1$days_incubated)
names_meanse <- rbind(paste0("mean day ", names_days), paste0("se day ", names_days))
colnames(my_mean_se2) <- names_meanse 

# Round the numbers to 2 decimal places and convert the numeric table to a data frame.
mean_se <- as.data.frame(round(my_mean_se2, digits = 2)); mean_se

#         mean day 0 se day 0 mean day 4 se day 4 mean day 10 se day 10 mean day 17
#i.x33.031      24.18    13.91       3.29     0.93        0.04      0.04       26.63
#i.x35.034       0.04     0.03       0.01     0.01        0.01      0.00        0.06
#i.x36.017       0.18     0.09       0.00     0.00        0.03      0.02        0.00
#i.x39.959       0.02     0.02       0.00     0.00        0.01      0.00        0.00
#i.x40.023       0.18     0.10       0.17     0.09        0.02      0.02        0.16
#          se day 17 mean day 24 se day 24 mean day 66 se day 66 mean day 81 se day 81
#i.x33.031     19.46      113.34     68.58       45.14     17.49       10.14      4.67
#i.x35.034      0.04        0.20      0.12        0.05      0.01        0.02      0.01
#i.x36.017      0.00        0.01      0.01        0.00      0.00        0.00      0.00
#i.x39.959      0.00        0.01      0.01        0.00      0.00        0.00      0.00
#i.x40.023      0.12        0.38      0.21        0.29      0.24        0.01      0.00
#          mean day 94 se day 94 mean day 116 se day 116
#i.x33.031       13.04      6.00         0.09       0.07
#i.x35.034        0.01      0.01         0.00       0.00
#i.x36.017        0.00      0.00         0.14       0.07
#i.x39.959        0.01      0.00         0.02       0.01
#i.x40.023        0.02      0.01         0.02       0.01

# Convert the data frame to a tibble with rowname set as a new column
as_tibble(mean_se, rownames = "rowname")

# A tibble: 5 x 19
#  rowname `mean day 0` `se day 0` `mean day 4` `se day 4` `mean day 10` `se day 10` `mean day 17` `se day 17` `mean day 24` `se day 24`
#  <chr>          <dbl>      <dbl>        <dbl>      <dbl>         <dbl>       <dbl>         <dbl>       <dbl>         <dbl>       <dbl>
#1 i.x33.~        24.2       13.9          3.29       0.93          0.04        0.04         26.6        19.5         113.         68.6 
#2 i.x35.~         0.04       0.03         0.01       0.01          0.01        0             0.06        0.04          0.2         0.12
#3 i.x36.~         0.18       0.09         0          0             0.03        0.02          0           0             0.01        0.01
#4 i.x39.~         0.02       0.02         0          0             0.01        0             0           0             0.01        0.01
#5 i.x40.~         0.18       0.1          0.17       0.09          0.02        0.02          0.16        0.12          0.38        0.21
# ... with 8 more variables: `mean day 66` <dbl>, `se day 66` <dbl>, `mean day 81` <dbl>, `se day 81` <dbl>, `mean day 94` <dbl>, `se day
#   94` <dbl>, `mean day 116` <dbl>, `se day 116` <dbl>

【讨论】：

我认为这些列的标签不正确。我认为这是 OP 遇到的同样问题，它们被排序为字符而不是数字。

【解决方案3】：

我认为你开始的方向是正确的。接下来，您可以获取长格式数据，按照您想要的方式更改列名，然后以宽格式获取数据。

library(dplyr)
library(tidyr)

df1 %>%
  group_by(days_incubated = as.numeric(days_incubated)) %>%
  summarise(across(i.x33.031:i.x40.023,  list(mean = mean, se = plotrix::std.error))) %>%
  pivot_longer(cols = -days_incubated, 
               names_to = c('mass', 'col'), 
               names_sep = '_') %>%
  #If you need column exactly as shown
  mutate(col = paste(col, 'day')) %>%
  pivot_wider(names_from = c(col, days_incubated), values_from = value, names_sep = ' ')


#  mass  `mean day 0` `se day 0` `mean day 4` `se day 4` `mean day 10` `se day 10`
#  <chr>        <dbl>      <dbl>        <dbl>      <dbl>         <dbl>       <dbl>
#1 i.x3…       24.2      13.9          45.1      17.5            3.29      0.932  
#2 i.x3…        0.036     0.0287        0.048     0.0116         0.014     0.00748
#3 i.x3…        0.178     0.0867        0         0              0         0      
#4 i.x3…        0.02      0.0155        0         0              0         0      
#5 i.x4…        0.182     0.0978        0.294     0.240          0.174     0.0941 
# … with 12 more variables: `mean day 17` <dbl>, `se day 17` <dbl>, `mean day
#   24` <dbl>, `se day 24` <dbl>, `mean day 66` <dbl>, `se day 66` <dbl>, `mean day
#   81` <dbl>, `se day 81` <dbl>, `mean day 94` <dbl>, `se day 94` <dbl>, `mean day
#   116` <dbl>, `se day 116` <dbl>

【讨论】：