Оценить экспоненциальное ограничение в степенном распределении - PullRequest
Новая
       0

Оценить экспоненциальное ограничение в степенном распределении

9 голосов
/ 30 января 2012

Поскольку я занимался анализом социальных сетей, я столкнулся с проблемой подбора распределения вероятностей по степени сети.

Итак, у меня есть распределение вероятностей P(X >= x), которое при визуальном осмотре следует степенному закону с экспоненциальным отсечением, а не чисто степенному закону (прямая линия).

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

f (x) = x ** альфа * exp (бета * x)

Как я могу оценить параметры alpha и beta, используя Python?

Я знаю, что пакет scipy.stats.powerlaw существует, и у него есть функция .fit(), но она, похоже, не справляется с работой, поскольку возвращает только местоположение и масштаб графика, что, по-видимому, полезно только для обычного распределение? В этом пакете также недостаточно учебников.

P.S. Мне хорошо известно о реализации CLauset и др. , но, похоже, они не обеспечивают способы оценки параметров альтернативных распределений.

Ответы [ 4 ]

3 голосов
/ 23 апреля 2012

Функция scipy.stats.powerlaw.fit все еще может работать для ваших целей. Это немного сбивает с толку, как работают дистрибутивы в scipy.stats (документация для каждого из них ссылается на необязательные параметры loc и scale, хотя не все из них используют эти параметры, и каждый использует их по-своему). Если вы посмотрите на документы:

http://docs.scipy.org/doc/scipy/reference/generated/scipy.stats.powerlaw.html

есть также второй необязательный параметр "a", который является "параметрами формы". В случае powerlaw, это содержит единственный параметр. Не беспокойтесь о "loc" и "scale".

Edit: Извините, забыл, что вы тоже хотели бета-параметр. Лучше всего вы можете определить нужную вам функцию powerlaw, а затем использовать общие алгоритмы подбора scipy для изучения параметров. Например: http://www.scipy.org/Cookbook/FittingData#head-5eba0779a34c07f5a596bbcf99dbc7886eac18e5

1 голос
/ 21 августа 2017

Powerlaw * Библиотека 1002 * может напрямую использоваться для оценки параметров следующим образом:

  1. Установить все зависимости питонов: 

    pip install powerlaw mpmath scipy

  2. Запуск подгонки пакета powerlaw в среде python: 

    import powerlaw
data = [5, 4, ... ]
results = powerlaw.Fit(data)

  3. получить параметры из результатов 

    results.truncated_power_law.parameter1 # power law  parameter (alpha)
results.truncated_power_law.parameter2 # exponential cut-off parameter (beta)
1 голос
/ 25 июня 2017

Вот средство оценки степени масштабирования и экспоненциальной скорости степенного закона с экспоненциальным отсечением путем максимизации вероятности в R: 

# Input: Data vector, lower threshold
# Output: List, giving type ("powerexp"), scaling exponent, exponential rate, lower threshold, log-likelihood


powerexp.fit <- function(data,threshold=1,method="constrOptim",initial_rate=-1) {
  x <- data[data>=threshold]
  negloglike <- function(theta) {
    -powerexp.loglike(x,threshold,exponent=theta[1],rate=theta[2])
  }
  # Fit a pure power-law distribution
  pure_powerlaw <- pareto.fit(data,threshold)
  # Use this as a first guess at the exponent
  initial_exponent <- pure_powerlaw$exponent
  if (initial_rate < 0) { initial_rate <- exp.fit(data,threshold)$rate }
  minute_rate <- 1e-6
  theta_0 <- as.vector(c(initial_exponent,initial_rate))
  theta_1 <- as.vector(c(initial_exponent,minute_rate))
  switch(method,
    constrOptim = {
      # Impose the constraint that rate >= 0
      # and that exponent >= -1
      ui <- rbind(c(1,0),c(0,1))
      ci <- c(-1,0)
      # Can't start with values on the boundary of the feasible set so add
      # tiny amounts just in case
      if (theta_0[1] == -1) {theta_0[1] <- theta_0[1] + minute_rate}
      if (theta_0[2] == 0) {theta_0[2] <- theta_0[2] + minute_rate}
      est <- constrOptim(theta=theta_0,f=negloglike,grad=NULL,ui=ui,ci=ci)
      alpha <- est$par[1]
      lambda <- est$par[2]
      loglike <- -est$value},
    optim = {
      est <- optim(par=theta_0,fn=negloglike)
      alpha <- est$par[1]
      lambda <- est$par[2]
      loglike <- -est$value},
    nlm = {
      est.0 <- nlm(f=negloglike,p=theta_0)
      est.1 <- nlm(f=negloglike,p=theta_1)
      est <- est.0
      if (-est.1$minimum > -est.0$minimum) { est <- est.1;cat("NLM had to switch\n") }
      alpha <- est$estimate[1]
      lambda <- est$estimate[2]
      loglike <- -est$minimum},
    {cat("Unknown method",method,"\n"); alpha<-NA; lambda<-NA; loglike<-NA}
  )
  fit <- list(type="powerexp", exponent=alpha, rate=lambda, xmin=threshold,
              loglike=loglike, samples.over.threshold=length(x))
  return(fit)
}

Проверьте https://github.com/jeffalstott/powerlaw/ для получения дополнительной информации

0 голосов
/ 22 сентября 2018

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

Дистрибутив, который я хотел приспособить, представляет собой закон обрезания степеней со сдвигом, такой как описанный в этой статье .Используя ту же запись, моя подгонка - 

f(x) = (x + x0)**alpha * exp(-beta*x)

, поэтому просто добавьте третий параметр x0 в ваш дистрибутив.Обратите внимание, что я предполагаю, что beta является положительным, и я просто вынимаю знак снаружи (я думаю, это проясняет, что ваша экспонента уменьшается).

Реализация выглядит следующим образом: 

import numpy as np.
import scipy.optimize as opt

def distribution(x, alpha, beta, x0):
    return (x + x0)**alpha * np.exp(-beta *x)

# ... I prepare my data here

fit = opt.curve_fit(distribution, x_data, y_data) # you can pass guess for the parameters/errors
alpha, beta, x0 = fit[0]

Это результат:

fit

...