Журнал Журнал Вероятностная Диаграмма в R
Я уверен, что это легко, но я рвал на себе волосы, пытаясь выяснить, как это сделать в R.
У меня есть некоторые данные, которые я пытаюсь приспособить к распределению степенного закона. Для этого вам необходимо вывести данные на кумулятивную диаграмму вероятности log-log. Ось Y - это журнал частоты данных (или, если хотите, логарифмическая вероятность), а ось X - это журнал значений. Если это прямая линия, то она соответствует распределению степенного закона, а градиент определяет параметр степенного закона.
Если мне нужна частота данных, я могу просто использовать функцию ecdf():
Мой набор данных называется Profits.negative, и это просто длинный список торговых прибылей, которые были меньше нуля (и я намеренно преобразовал их все в положительные числа, чтобы избежать проблем с регистрацией в дальнейшем).
Так что я могу напечатать
plot(ecdf(Profits.negative))
И я получаю удобную эмпирическую функцию CDF. Все, что мне нужно сделать, это преобразовать обе оси в логарифмические шкалы. Я могу сделать ось X:
Profits.negative.logs <- log(Profits.negative)
plot(ecdf(Profits.negative.logs))
Почти готово! Мне просто нужно разобраться, как войти по оси Y! Но я не могу этого сделать и не могу понять, как извлечь фигуры из объекта ecdf. Кто-нибудь может помочь?
Я знаю, что есть функция power.law.fit, но она только оценивает параметры - я хочу построить данные и посмотреть, совпадают ли они.
3 ответа
Вы можете подобрать и построить степенные законы, используя пакет poweRlaw. Вот пример. Сначала мы генерируем некоторые данные из дистрибутива с тяжелыми хвостами:
set.seed(1)
x = round(rlnorm(100, 3, 2)+1)
Затем мы загружаем пакет и создаем объект данных и объект displ:
library(poweRlaw)
m = displ$new(x)
Мы можем оценить xmin
и параметр масштабирования:
est = estimate_xmin(m))
и установить параметры
m$setXmin(est[[2]])
m$setPars(est[[3]])
Затем нанесите данные и добавьте подходящую линию:
plot(m)
lines(m, col=2)
Получить:
Генерация данных в первую очередь (насамом деле вы расстались;)):
set.seed(1)
Profits.negative <- runif(1e3, 50, 100) + rnorm(1e2, 5, 5)
Ведение журнала и ecdf
:
Profits.negative.logs <- log(Profits.negative)
fn <- ecdf(Profits.negative.logs)
ecdf
возвращает функцию, и если вы хотите извлечь что-то из нее - это хорошая идея, чтобы посмотреть на закрытие функции:
ls(environment(fn))
# [1] "f" "method" "n" "nobs" "x" "y" "yleft" "yright"
Ну, теперь мы можем получить доступ x
а также y
:
x <- environment(fn)$x
y <- environment(fn)$y
Вероятно, это то, что вам нужно. В самом деле, plot(fn)
а также plot(x,y,type="l")
показать практически те же результаты. Для входа по оси Y вам нужно всего лишь:
plot(x,log(y),type="l")
Вот подход с использованием ggplot2
:
library(ggplot2)
# data
set.seed(1)
x = round(rlnorm(100, 3, 2)+1)
# organize data into a df
df <- data.frame(x = sort(x, decreasing = T),
pk <- ecdf(x)(x),
k <- seq_along(x))
# plot
ggplot(df, aes(x=k, y= pk)) + geom_point(alpha=0.5) +
coord_trans(x = 'log10', y = 'log10') +
scale_x_continuous(breaks = trans_breaks("log10", function(x) 10^x), labels = trans_format("log10", math_format(10^.x))) +
scale_y_continuous(breaks = trans_breaks("log10", function(x) 10^x), labels = trans_format("log10", math_format(10^.x)))