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

Question

Я ищу что-то похожее по форме на weighted.mean(). Я нашел некоторые решения с помощью поиска, которые выписывают всю функцию, но были бы признательны за что-то более удобное для пользователя.

Answer 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)

Answer 2

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

median(rep(x, times=w))

Answer 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() является надежным решением

Answer 4

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

Answer 5

Публикация исходного кода для функций spatstat (упомянутых в ответе пользователя 2522202) здесь, потому что люди могут не захотеть устанавливать этот пакет, который имеет много зависимостей, просто для получения взвешенного медианы / квантилей , Сами функции не имеют зависимостей. Я добавил код 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))
}

Answer 6

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

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[w>0] ~ x[w>0], 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[w>0] ~ x[w>0], 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[w>0] ~ x[w>0], 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

Answer 7

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

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

Answer 8

Можно также использовать 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.