Надежно получить обратную функцию квантиля

Question

Я читал другие посты (например, здесь ) о получении "обратного" квантиля - то есть о том, чтобы получить процентиль, соответствующий определенному значению в серии значений.

Однако ответы не дают мне того же значения, что и квантиль для того же ряда данных .

Я также исследовал, что квантиль предоставляет 9 различных алгоритмов для вычисления процентиля.

Итак, мой вопрос: есть ли надежный способ получить обратную функцию квантиля? ecdf не принимает аргумент типа, так что не похоже, чтобы они использовали один и тот же метод.

Воспроизводимый пример:

# Simple data
x = 0:10
pcntile = 0.5


# Get value corresponding to a percentile using quantile
(pcntile_value <- quantile(x, pcntile))     

# 50%    
# 5               # returns 5 as expected for 50% percentile     



# Get percentile corresponding to a value using ecdf function
(pcntile_rev <- ecdf(x)(5))                


# [1] 0.5454545   #returns 54.54% as the percentile for the value 5


# Not the same answer as quantile produces

Answer 1

Ответ в ссылке действительно хорош, но, возможно, это поможет, взглянуть на ecdf Просто запустите следующий код:

# Simple data
x = 0:10
p0 = 0.5

# Get value corresponding to a percentile using quantile
sapply(c(1:7), function(i) quantile(x, p0, type = i))
# 50% 50% 50% 50% 50% 50% 50% 
# 5.0 5.0 5.0 4.5 5.0 5.0 5.0 

Таким образом, это не вопрос типа. Вы можете войти в функцию, используя debug:

# Get percentile corresponding to a value using ecdf function
debug(ecdf)
my_ecdf <- ecdf(x)

Важнейшей частью является

rval <- approxfun(vals, cumsum(tabulate(match(x, vals)))/n, 
    method = "constant", yleft = 0, yright = 1, f = 0, ties = "ordered")

После этого вы можете проверить

data.frame(x = vals, y = round(cumsum(tabulate(match(x, vals)))/n, 3), stringsAsFactors = FALSE)

и если вы поделитесь на n=11, результат неудивителен Как уже говорилось, для теории взгляните на другой ответ.

Кстати, вы также можете построить график функции

plot(my_ecdf)

По поводу вашего комментария. Я думаю, что это не вопрос надежности, а вопрос о том, как определить «обратную функцию распределения, если она не существует»:

Хороший справочник по обобщенным инверсиям: Пол Эмбрехтс, Мариус Хоферт: «Замечание по обобщенным инверсиям», Math Meth Oper Res (2013) 77: 423–432 DOI

Answer 2

ecdf дает результат формулы в документации.

x <- 0:10
Fn <- ecdf(x)

Теперь объект Fn является интерполяционной пошаговой функцией.

str(Fn)
#function (v)  
# - attr(*, "class")= chr [1:3] "ecdf" "stepfun" "function"
# - attr(*, "call")= language ecdf(x)

И он сохраняет исходные x значения и соответствующие y значения.

environment(Fn)$x
# [1]  0  1  2  3  4  5  6  7  8  9 10

environment(Fn)$y
# [1] 0.09090909 0.18181818 0.27272727 0.36363636 0.45454545 0.54545455
# [7] 0.63636364 0.72727273 0.81818182 0.90909091 1.00000000

Последние являются точно такими же значениями, поскольку в результате в документации указывается формула, используемая для их вычисления. От help('ecdf'):

Для наблюдений x = (x1, x2, ... xn), Fn - это доля
наблюдения меньше или равны t, т. е.

Fn (t) = # {xi <= t} / n = 1 / n сумма (i = 1, n) Индикатор (xi <= t). </p>

Вместо 1:length(x) я буду использовать seq_along.

seq_along(x)/length(x)
# [1] 0.09090909 0.18181818 0.27272727 0.36363636 0.45454545 0.54545455
# [7] 0.63636364 0.72727273 0.81818182 0.90909091 1.00000000
Fn(x)
# [1] 0.09090909 0.18181818 0.27272727 0.36363636 0.45454545 0.54545455
# [7] 0.63636364 0.72727273 0.81818182 0.90909091 1.00000000