plfit/plfit_v1.py
def plfit(x,nosmall=False,finite=False):
    """
    A Python implementation of the Matlab code http://www.santafe.edu/~aaronc/powerlaws/plfit.m
    from http://www.santafe.edu/~aaronc/powerlaws/

    See A. Clauset, C.R. Shalizi, and M.E.J. Newman, "Power-law distributions
    in empirical data" SIAM Review, to appear (2009). (arXiv:0706.1062)
    http://arxiv.org/abs/0706.1062
    """
    xmins = unique(x)
    xmins = xmins[1:-1]
    dat = xmins * 0 
    z = sort(x)
    for xm in arange(len(xmins)):
        xmin = xmins[xm]
        z    = z[z>=xmin] 
        n    = float(len(z))
        # estimate alpha using direct MLE
        a    =  n / sum( log(z/xmin) )
        if nosmall:
            # 4. For continuous data, PLFIT can return erroneously large estimates of 
            #    alpha when xmin is so large that the number of obs x >= xmin is very 
            #    small. To prevent this, we can truncate the search over xmin values 
            #    before the finite-size bias becomes significant by calling PLFIT as
            if (a-1)/sqrt(n) > 0.1:
                #dat(xm:end) = [];
                dat = dat[:xm]
                xm = len(xmins)+1
                break
        # compute KS statistic
        cx   = arange(n)/float(n)  #data
        cf   = 1-(xmin/z)**a  # fitted
        dat[xm] = max( abs(cf-cx) )
    D     = min(dat);
    #xmin  = xmins(find(dat<=D,1,'first'));
    xmin  = xmins[argmin(dat)]
    z     = x[x>=xmin]
    n     = len(z)
    alpha = 1 + n / sum( log(z/xmin) )
    if finite:
        alpha = alpha*(n-1)/n+1/n
    if n < 50 and ~finite:
        print '(PLFIT) Warning: finite-size bias may be present.'
    L = n*log((alpha-1)/xmin) - alpha*sum(log(z/xmin));
    return xmin,alpha,L,dat
def plfit(x,nosmall=False,finite=False):
    """
    A Python implementation of the Matlab code http://www.santafe.edu/~aaronc/powerlaws/plfit.m
    from http://www.santafe.edu/~aaronc/powerlaws/

    See A. Clauset, C.R. Shalizi, and M.E.J. Newman, "Power-law distributions
    in empirical data" SIAM Review, to appear (2009). (arXiv:0706.1062)
    http://arxiv.org/abs/0706.1062
    """
    xmins = unique(x)
    xmins = xmins[1:-1]
    dat = xmins * 0 
    z = sort(x)
    for xm in arange(len(xmins)):
        xmin = xmins[xm]
        z    = z[z>=xmin] 
        n    = float(len(z))
        # estimate alpha using direct MLE
        a    =  n / sum( log(z/xmin) )
        if nosmall:
            # 4. For continuous data, PLFIT can return erroneously large estimates of 
            #    alpha when xmin is so large that the number of obs x >= xmin is very 
            #    small. To prevent this, we can truncate the search over xmin values 
            #    before the finite-size bias becomes significant by calling PLFIT as
            if (a-1)/sqrt(n) > 0.1:
                #dat(xm:end) = [];
                dat = dat[:xm]
                xm = len(xmins)+1
                break
        # compute KS statistic
        cx   = arange(n)/float(n)  #data
        cf   = 1-(xmin/z)**a  # fitted
        dat[xm] = max( abs(cf-cx) )
    D     = min(dat);
    #xmin  = xmins(find(dat<=D,1,'first'));
    xmin  = xmins[argmin(dat)]
    z     = x[x>=xmin]
    n     = len(z)
    alpha = 1 + n / sum( log(z/xmin) )
    if finite:
        alpha = alpha*(n-1)/n+1/n
    if n < 50 and ~finite:
        print '(PLFIT) Warning: finite-size bias may be present.'
    L = n*log((alpha-1)/xmin) - alpha*sum(log(z/xmin));
    return xmin,alpha,L,dat
A Python implementation of the Matlab code http://www.santafe.edu/~aaronc/powerlaws/plfit.m
    from http://www.santafe.edu/~aaronc/powerlaws/

    See A. Clauset, C.R. Shalizi, and M.E.J. Newman, "Power-law distributions
    in empirical data" SIAM Review, to appear (2009). (arXiv:0706.1062)
    http://arxiv.org/abs/0706.1062