Есть ли функция weighted.median()?

Ответ 1

В следующих пакетах есть функция для вычисления взвешенной медианы: "aroma.light", "isotone", "limma", "cwhmisc", "ergm", "laeken", "matrixStats", "PSCBS" и "bigvis" (на github).

Чтобы найти их, я использовал неоценимый findFn() в пакете sos, который является расширением для встроенной справки R.

findFn('weighted median')

Или

???'weighted median'

как??? является ярлыком таким же образом ?some.function для help(some.function)

Ответ 2

Чтобы вычислить взвешенную медиану вектора x, используя вектор одинаковой длины (целых) весов w:

median(rep(x, times=w))

Ответ 3

Некоторые примеры использования ответов от @wkmor1 и @Jaitropmange.

Я проверил 3 функции из 3 пакетов, isotone, laeken и matrixStats. Только matrixStats работает правильно. Другие два (как решение median(rep(x, times=w)) дают целочисленный вывод. Пока я вычислял средний возраст популяций, значение десятичных знаков имеет значение.

Воспроизводимый пример. Расчет среднего возраста популяции

df <- data.frame(age = 0:100,
                 pop = spline(c(4,7,9,8,7,6,4,3,2,1),n = 101)$y)

library(isotone)
library(laeken)
library(matrixStats)

isotone::weighted.median(df$age,df$pop)
# [1] 36
laeken::weightedMedian(df$age,df$pop)
# [1] 36
matrixStats::weightedMedian(df$age,df$pop)
# [1] 36.164
median(rep(df$age, times=df$pop))
# [1] 35

Резюме

matrixStats::weightedMedian() является надежным решением

Ответ 4

Действительно старый пост, но я только что наткнулся на него и провел несколько испытаний различных методов. spatstat::weighted.median(), по-видимому, примерно в 14 раз быстрее, чем median(rep(x, times=w)) и это действительно заметно, если вы хотите запустить функцию более пары раз. Тестирование было проведено с относительно большим опросом, около 15 000 человек.

Ответ 5

Если вы работаете с пакетом survey, предполагая, что вы определили его структуру опроса, а x представляет вашу переменную интереса:

svyquantile(~x, mydesign, c(0.5))

Ответ 6

Публикация исходного кода для функций spatstat (упомянутое в ответе user2522202) здесь, потому что люди могут не захотеть устанавливать этот пакет, который имеет много зависимостей, просто чтобы получить взвешенное медиану/квантили. Сами функции не имеют зависимостей. Я добавил код Roxygen на тот случай, если вы хотите поместить его в пакет.

#' Weighted quantile
#'
#' Function copied from **spatstat** package.
#'
#' @param x Vector of values
#' @param w Vector of weights
#' @param probs Vector of probabilities
#' @param na.rm Ignore missing data?
#' @export
weighted.quantile <- function(x, w, probs=seq(0,1,0.25), na.rm=TRUE) {
  x <- as.numeric(as.vector(x))
  w <- as.numeric(as.vector(w))
  if(anyNA(x) || anyNA(w)) {
    ok <- !(is.na(x) | is.na(w))
    x <- x[ok]
    w <- w[ok]
  }
  stopifnot(all(w >= 0))
  if(all(w == 0)) stop("All weights are zero", call.=FALSE)
  #'
  oo <- order(x)
  x <- x[oo]
  w <- w[oo]
  Fx <- cumsum(w)/sum(w)
  #'
  result <- numeric(length(probs))
  for(i in seq_along(result)) {
    p <- probs[i]
    lefties <- which(Fx <= p)
    if(length(lefties) == 0) {
      result[i] <- x[1]
    } else {
      left <- max(lefties)
      result[i] <- x[left]
      if(Fx[left] < p && left < length(x)) {
        right <- left+1
        y <- x[left] + (x[right]-x[left]) * (p-Fx[left])/(Fx[right]-Fx[left])
        if(is.finite(y)) result[i] <- y
      }
    }
  }
  names(result) <- paste0(format(100 * probs, trim = TRUE), "%")
  return(result)
}


#' Weighted median
#'
#' Function copied from **spatstat** package.
#'
#' @param x Vector of values
#' @param w Vector of weights
#' @param na.rm Ignore missing data?
#' @export
weighted.median <- function(x, w, na.rm=TRUE) {
  unname(weighted.quantile(x, probs=0.5, w=w, na.rm=na.rm))
}

Ответ 7

Можно также использовать stats::density для создания взвешенного PDF, а затем преобразовать его в CDF, как описано здесь:

my_wtd_q = function(x, w, prob, n = 4096) 
  with(density(x, weights = w/sum(w), n = n), 
       x[which.max(cumsum(y*(x[2L] - x[1L])) >= prob)])

Тогда my_wtd_q(x, w,.5) будет взвешенной медианой.

Можно также быть более осторожным, чтобы обеспечить повторную нормализацию общей площади под density.

Ответ 8

Используя источник из Deleet и данные из Икашницкого, средневзвешенную медиану можно рассчитать на основе:

df <- data.frame(age = 0:100,
                 pop = spline(c(4,7,9,8,7,6,4,3,2,1),n = 101)$y)

medianWeighted <- function(x, w) {
  x <- aggregate(w ~ x, FUN=sum)
  approxfun(filter(c(0,cumsum(x$w)/sum(x$w)), c(.5,.5), sides=1)[-1], x$x)(.5)
}
medianWeighted(df$age,df$pop) #Interpolates between observed Numbers
#[1] 36.164

medianWeightedI <- function(x, w) { 
  w <- w[order(x)]
  x <- x[order(x)]
  x[which.min(abs(filter(c(0,cumsum(w)/sum(w)), c(.5,.5), sides=1)[-1] - 0.5))]
}
medianWeightedI(df$age,df$pop) #Takes only numbers which have been observed
#[1] 36

В случае, если вы также хотите рассчитать взвешенные квантили.

quantileWeighted <- function(x, w, probs = seq(0, 1, 0.25)) {
  x <- aggregate(w ~ x, FUN=sum)
  approxfun(filter(c(0,cumsum(x$w)/sum(x$w)), c(.5,.5), sides=1)[-1], x$x, rule=2)(probs)
}
quantileWeighted(df$age, df$pop)
#[1]   0.00000  20.21336  36.16400  55.98371 100.00000

quantileWeightedI <- function(x, w, probs = seq(0, 1, 0.25)) {
  x <- aggregate(w ~ x, FUN=sum)
  stepfun(cumsum(x$w[-nrow(x)])/sum(x$w[-nrow(x)]), x$x)(probs)
}
quantileWeightedI(df$age, df$pop)
#[1]   0  20  36  56 100