几何平均值：有内置的吗？答案

【问题标题】：Geometric Mean: is there a built-in?几何平均值：有内置的吗？
【发布时间】：2011-02-05 20:26:21
【问题描述】：

我试图找到几何平均值的内置函数，但找不到。

（显然，在 shell 中工作时，内置不会为我节省任何时间，我也不怀疑准确性有任何差异；对于脚本，我尝试尽可能多地使用内置，其中（累积）性能提升通常很明显。

如果没有（我怀疑是这种情况），这是我的。

gm_mean = function(a){prod(a)^(1/length(a))}

【问题讨论】：

小心负数和溢出。 prod(a) 将很快低于或溢出。我尝试使用一个大列表来计时，并使用您的方法快速获得 Inf 与 exp(mean(log(x))); 的 1.4；舍入问题可能非常严重。
我只是快速编写了上面的函数，因为我确信在发布这个 Q 5 分钟后，有人会告诉我 R 的 gm 内置函数。所以没有内置的，所以根据你的评论花时间重新编码是肯定的。 + 1 来自我。
我刚刚标记了这个 geometric-mean 和 built-in，9 年后。

标签： r statistics built-in geometric-mean

【解决方案1】：

没有，但是有少数人写过，比如here。

另一种可能性是使用这个：

exp(mean(log(x)))

【讨论】：

使用 exp(mean(log(x))) 的另一个优点是您可以处理很长的大数列表，这在使用 prod() 使用更明显的公式时会出现问题。请注意 prod(a)^(1/length(a)) 和 exp(mean(log(a))) 给出相同的答案。
链接已修复

【解决方案2】：

我完全按照马克所说的来使用。这样，即使使用tapply，您也可以使用内置的mean 功能，无需定义您的！例如，计算 data$value 的每组几何平均值：

exp(tapply(log(data$value), data$group, mean))

【讨论】：

【解决方案3】：

我们可以使用psych package并调用geometric.mean函数。

【讨论】：

psych::geometric.mean()
我会说，这些功能应该是系列而不是它们的增长，至少作为一种选择。

【解决方案4】：

exp(mean(log(x)))

除非 x 中有 0，否则将起作用。如果是这样，日志将产生 -Inf (-Infinite)，它总是导致几何平均值为 0。

一种解决方案是在计算平均值之前删除 -Inf 值：

geo_mean <- function(data) {
    log_data <- log(data)
    gm <- exp(mean(log_data[is.finite(log_data)]))
    return(gm)
}

您可以使用单行来执行此操作，但这意味着计算日志两次，效率低下。

exp(mean(log(i[is.finite(log(i))])))

【讨论】：

为什么可以计算两次日志：exp(mean(x[x!=0]))
这两种方法的均值都是错误的，因为如果您过滤 x 然后将其传递给 mean，则均值的分母 sum(x) / length(x) 是错误的。
我认为过滤是个坏主意，除非您明确表示要这样做（例如，如果我正在编写通用函数，我不会将过滤设为默认值）-好的，如果这是一次性代码，并且您已经非常仔细地考虑过在您的问题上下文中过滤清零的实际含义（！）
根据定义，一组包含零的数字的几何平均值应该为零！ math.stackexchange.com/a/91445/221143

【解决方案5】：

如果您的数据中存在缺失值，这种情况并不少见。您需要再添加一个参数。

你可以试试下面的代码：

exp(mean(log(i[ is.finite(log(i)) ]), na.rm = TRUE))

【讨论】：

【解决方案6】：

这是一个矢量化、零和 NA 容错函数，用于计算 R 中的几何平均值。对于 x 包含非正值的情况，涉及 length(x) 的详细 mean 计算是必要的。

gm_mean = function(x, na.rm=TRUE){
  exp(sum(log(x[x > 0]), na.rm=na.rm) / length(x))
}

感谢@ben-bolker 注意到na.rm 传递和@Gregor 确保它正常工作。

我认为某些 cmets 与数据和零中 NA 值的错误等效有关。在我想到的应用程序中，它们是相同的，但当然这通常不是真的。因此，如果您想包含可选的零传播，并在删除 NA 的情况下区别对待 length(x)，以下是上述函数的稍长替代方案。

gm_mean = function(x, na.rm=TRUE, zero.propagate = FALSE){
  if(any(x < 0, na.rm = TRUE)){
    return(NaN)
  }
  if(zero.propagate){
    if(any(x == 0, na.rm = TRUE)){
      return(0)
    }
    exp(mean(log(x), na.rm = na.rm))
  } else {
    exp(sum(log(x[x > 0]), na.rm=na.rm) / length(x))
  }
}

请注意，它还会检查是否有任何负值，并返回一个信息更丰富、更合适的NaN，说明几何平均值不是为负值定义的（而是为零定义的）。感谢一直关注我的情况的评论者。

【讨论】：

将na.rm 作为参数传递不是更好吗（即让用户决定他们是否要容忍NA，以与其他R摘要函数保持一致）？我对自动排除零感到紧张——我也会选择这样做。
也许您将na.rm 作为选项传递是对的。我会更新我的答案。至于排除零，几何平均值对于非正值是未定义的，包括零。以上是几何平均值的常见修复，其中零（或在这种情况下所有非零）被赋予虚拟值 1，这对乘积没有影响（或等效地，对数和中的零）。
您的na.rm 传递无法按编码工作...请参阅gm_mean(c(1:3, NA), na.rm = T)。需要从向量子集中去掉& !is.na(x)，由于sum的第一个arg是...，所以需要按名称传递na.rm = na.rm，还需要排除0的和@ 987654340@ 来自 length 调用中的向量。
当心：x 仅包含零，例如 x <- 0，exp(sum(log(x[x>0]), na.rm = TRUE)/length(x)) 将 1 用于几何平均值，这是没有意义的。
假设 na.rm = TRUE，不是必须是 length(x[!is.na(x) & x > 0]) 吗？

【解决方案7】：

EnvStats package 具有 geoMean 和 geoSd 的功能。

【讨论】：

【解决方案8】：

这个版本提供了比其他答案更多的选项。

它允许用户区分不是（实）数字的结果和不可用的结果。如果存在负数，则答案将不是实数，因此返回 NaN。如果都是NA 值，那么该函数将返回NA_real_，以反映实际值实际上不可用。这是一个微妙的差异，但可能会产生（稍微）更可靠的结果。
第一个可选参数zero.rm 旨在允许用户让零影响输出而不使其为零。如果zero.rm 设置为FALSE 并且eta 设置为NA_real_（其默认值），则零具有将结果缩小为一的效果。我对此没有任何理论依据 - 似乎更有意义的是不忽略零，而是“做一些不涉及自动使结果为零的事情”。
eta 是一种处理零的方法，其灵感来自以下讨论：https://support.bioconductor.org/p/64014/

geomean <- function(x,
                    zero.rm = TRUE,
                    na.rm = TRUE,
                    nan.rm = TRUE,
                    eta = NA_real_) {
    nan.count <- sum(is.nan(x))
     na.count <- sum(is.na(x))
  value.count <- if(zero.rm) sum(x[!is.na(x)] > 0) else sum(!is.na(x))

  #Handle cases when there are negative values, all values are missing, or
  #missing values are not tolerated.
  if ((nan.count > 0 & !nan.rm) | any(x < 0, na.rm = TRUE)) {
    return(NaN)
  }
  if ((na.count > 0 & !na.rm) | value.count == 0) {
    return(NA_real_)
  }

  #Handle cases when non-missing values are either all positive or all zero.
  #In these cases the eta parameter is irrelevant and therefore ignored.
  if (all(x > 0, na.rm = TRUE)) {
    return(exp(mean(log(x), na.rm = TRUE)))
  }
  if (all(x == 0, na.rm = TRUE)) {
    return(0)
  }

  #All remaining cases are cases when there are a mix of positive and zero
  #values.
  #By default, we do not use an artificial constant or propagate zeros.
  if (is.na(eta)) {
    return(exp(sum(log(x[x > 0]), na.rm = TRUE) / value.count))
  }
  if (eta > 0) {
    return(exp(mean(log(x + eta), na.rm = TRUE)) - eta)
  }
  return(0) #only propagate zeroes when eta is set to 0 (or less than 0)
}

【讨论】：

您能否添加一些详细信息来解释这与现有解决方案有何不同/如何改进？（除非必要，否则我个人不想为这样的实用程序添加像 dplyr 这样的重度依赖项......）
我同意，case_whens 有点傻，所以我删除了它们和对ifs 的依赖。我还提供了一些详细说明。
我同意你的后一个想法，将nan.rm 的默认值更改为TRUE 以对齐所有三个``.rm` 参数。
另一种风格挑剔。 ifelse 专为矢量化而设计。如果要检查一个条件，使用value.count <- if(zero.rm) sum(x[!is.na(x)] > 0) else sum(!is.na(x)) 会更惯用
它看起来也比ifelse 好。改变了。谢谢！

【解决方案9】：

exp(mean(log(x1))) == prod(x1)^(1/length(x1))

【讨论】：