У меня есть несколько частот запросов, и мне нужно оценить коэффициент закона Ципфа. Это верхние частоты:
26486
12053
5052
3033
2536
2391
1444
1220
1152
1039
У меня есть несколько частот запросов, и мне нужно оценить коэффициент закона Ципфа. Это верхние частоты:
26486
12053
5052
3033
2536
2391
1444
1220
1152
1039
Ответы:
Обновление Я обновил код с оценкой максимального правдоподобия согласно предложению @whuber. Минимизация суммы квадратов различий между логарифмическими теоретическими вероятностями и логарифмическими частотами, хотя дает ответ, была бы статистической процедурой, если бы можно было показать, что это своего рода М-оценка. К сожалению, я не мог придумать ни одного, который мог бы дать такие же результаты.
Вот моя попытка. Я вычисляю логарифмы частот и пытаюсь привести их в соответствие с логарифмами теоретических вероятностей, которые дает эта формула . Конечный результат кажется разумным. Вот мой код в R.
Тогда наилучшее квадратичное соответствие будет .s = 1,47
Максимальное правдоподобие в R можно выполнить с помощью
mle
функции (изstats4
пакета), которая услужливо рассчитывает стандартные ошибки (если указана правильная функция отрицательного максимального правдоподобия):Вот график соответствия в масштабе log-log (опять же, как предложил @whuber):
Красная линия - это сумма подходящих квадратов, зеленая линия - соответствие с максимальным правдоподобием.
источник
Перед любой проблемой оценки стоит несколько вопросов :
Оцените параметр.
Оцените качество этой оценки.
Изучите данные.
Оцените подгонку.
Для тех, кто использует статистические методы для понимания и общения, первое никогда не должно делаться без других.
Таким образом, логарифмическая вероятность для данных
Учитывая природу закона Ципфа, правильный способ составить график этого подбора - на графике log-log , где подгонка будет линейной (по определению):
Чтобы оценить правильность подбора и изучить данные, посмотрите на остатки (данные / подбор, снова ось log-log):
Поскольку остатки кажутся случайными, в некоторых приложениях мы могли бы принять закон Ципфа (и нашу оценку параметра) в качестве приемлемого, хотя и грубого описания частот . Этот анализ показывает, однако, что было бы ошибкой предполагать, что эта оценка имеет какое-либо объяснительное или прогнозное значение для набора данных, рассматриваемого здесь.
источник
Один из вероятностных языков программирования, такой как PyMC3, делает эту оценку относительно простой. Другие языки включают Stan, у которого есть отличные возможности и поддерживающее сообщество.
Вот моя реализация Python модели, адаптированной к данным OP (также на Github ):
Чтобы обеспечить некоторую базовую диагностику выборки, мы можем видеть, что выборка «хорошо перемешивалась», так как мы не видим никакой структуры в трассе:
Для запуска кода необходим Python с установленными пакетами Theano и PyMC3.
Спасибо @ w-huber за отличный ответ и комментарии!
источник
Вот моя попытка привести данные в соответствие, оценить и изучить результаты с помощью VGAM:
В нашем случае нулевая гипотеза Чи-квадрата состоит в том, что данные распределяются по закону Зипфа, поэтому большие значения р подтверждают утверждение, что данные распределяются в соответствии с ним. Обратите внимание, что даже очень большие значения p не являются доказательством, это просто показатель.
источник
Опять же, UWSE предоставляет только непротиворечивую оценку - без доверительных интервалов, и мы можем увидеть некоторый компромисс в точности. Приведенное выше решение mpiktas также является приложением UWSE - хотя программирование необходимо. Для полного объяснения оценки см .: https://paradsp.wordpress.com/ - полностью внизу.
источник
Мое решение пытается дополнить ответы, предоставленные mpiktas и whuber, которые осуществляют реализацию на Python. Наши частоты и диапазоны х:
Поскольку наша функция не определена во всем диапазоне, нам нужно проверять, что мы нормализуем каждый раз, когда вычисляем ее. В дискретном случае простым приближением является деление на сумму всех y (x). Таким образом, мы можем сравнить различные параметры.
Результат дает нам наклон 1.450408, как и в предыдущих ответах.
источник