**Power-law Distributions in Empirical Data**

This page is a companion for the *SIAM
Review* paper on power-law distributions in empirical data, written by
Aaron Clauset (me),
Cosma R. Shalizi and
M.E.J. Newman.
This page hosts implementations of the methods we describe in the article,
including several by authors other than us. Our goal is for the methods to
be widely accessible to the community. Python users may want to consider the
powerlaw package by Alstott et al.

*NOTE: we cannot provide technical support for code not written by us, and
we are busy with other projects now and so may not provide support for our own
code.*

**Journal References**

A. Clauset, C.R. Shalizi, and M.E.J. Newman, "Power-law distributions in empirical data" *SIAM Review* **51**(4), 661-703 (2009).
(arXiv:0706.1062, doi:10.1137/070710111)

Y. Virkar and A. Clauset, Power-law distributions in binned empirical data.
*Annals of Applied Statistics* **8**(1), 89 - 119 (2014). (arXiv:1208.3524; get the code)

**Random number generators**

This function generates continuous values randomly distributed according to
one of the five distributions discussed in the article (power law,
exponential, log-normal, stretched exponential, and power law with cutoff).
Usage information is included in the file; type 'help randht' at the Matlab
prompt for more information.

randht.m (Matlab, by Aaron Clauset)

randht.py (Python, by Joel Ornstein)

**Fitting a power-law distribution**

This function implements both the discrete and continuous maximum likelihood
estimators for fitting the power-law distribution to data, along with the
goodness-of-fit based approach to estimating the lower cutoff for the scaling
region. Usage information is included in the file; type 'help plfit' at the
Matlab prompt for more information.

plfit.m (Matlab, by Aaron Clauset)

plfit.r (R, by Laurent Dubroca)

plfit.py (Python, by Adam Ginsburg)

plfit.c (C++, by Wim Otte; includes plvar.c)

plfit.c (C++, by Tamas Nepusz)

plfit.py (Python, by Joel Ornstein)

**Visualizing the fitted distribution**

After several requests, I've written this function, which plots (on log-log
axes) the empirical distribution along with the fitted power-law distribution.
Usage information is included in the file; type 'help plplot' at the Matlab
prompt for more information.

plplot.m (Matlab, by Aaron Clauset)

plplot.py (Python, by Joel Ornstein)

**Estimating uncertainty in the fitted parameters**

This function implements the nonparametric approach for estimating the
uncertainty in the estimated parameters for the power-law fit found by the
plfit function. It too implements both continuous and discrete versions. Usage
information is included in the file; type 'help plvar' at the Matlab prompt
for more information.

plvar.m (Matlab, by Aaron Clauset)

plvar.c (C++, by Wim Otte; includes plfit.c)

plvar.py (Python, by Joel Ornstein)

**Calculating p-value for fitted power-law model**

This function implements the Kolmogorov-Smirnov test (which computes a

*p*-value for the estimated power-law fit to the data) for the power-law model. As above, it too implements both continuous and discrete versions of the test. Usage information is included in the file; type 'help plpva' at the Matlab prompt for more information.

plpva.m (Matlab, by Aaron Clauset)

plpva.r (R, by Laurent Dubroca; modified by Neal Walfield)

parplpva2.m (Matlab, by Casper Peterson, uses Parallel Toolbox)

plpva.py (Python, by Joel Ornstein)

**Riemann Zeta function**

The discrete estimator needs to calculate the Hurwitz Zeta function for
normalization. Matlab includes this function in the Symbolic Math Toolbox
(but be warned that their implementation becomes unstable for large alpha and
xmin, e.g., alpha>7 with xmin>150). There are also free versions available if
you don't have this toolbox. For instance, Paul Godfrey's
special functions library (via Matlab Central File Exchange) gives one, which we
mirror here (note, you need both these files; tip to Will Tracy).

deta.m (Matlab, by Paul Godfrey)

zeta.m (Matlab, by Paul Godfrey)

**Calculating likelihood-ratio test results**

The functions necessary to compute the log likelihood ratio tests is
implemented in the statistical programming
language R. Documentation of these functions is given in a separate file,
and the R functions themselves are in a downloadable tgz file (note: this is
not a proper R package, yet).

Documentation

R code (by Cosma Shalizi)

**Download all files**

Get the most up-to-date versions of the complete implementations.

Download all Matlab and R files (by Aaron Clauset and Cosma Shalizi)

Download Python package (by Jeff Alstott)

Download Python (2.6) package (by Javier del Molino Matamala)

Download Java package (by Peter Bloem)

Download R package (by Colin Gillespie)

**A note Matlab compatibility**

For the Matlab functions written by me (Aaron), all of them were designed to
be compatible with Matlab v7. They are *not* necessarily compatible with
older versions of Matlab. That being said, it should be possible to make them
compatible as the core functionality does not depend on v7 features.

**A note about bugs and alternative implementations**

The code provided here is provided as-is, with no warranty, with no guarantees
of technical support or maintenance, etc. If you experience problems while
using the code, please let the author(s) know via email. I am happy to host
(or link to) implementations of any of these functions in other programming
languages, in the interest of facilitating their more widespread use.
However, I cannot provide any technical support for that code.
The original pl* functions (Matlab) were written by Aaron Clauset and the LRT
functions (R) were written by Cosma Shalizi; all other language
implementations were written by members of the wider community.

Finally, if you use our code in an academic publication, it would be
courteous of you to thank me (Aaron) and Cosma in your acknowledgements for
providing you with implementations of the methods. If you use the
implementations by other authors, you should acknowledge them instead.

**A note the data sets**

The 24 data sets we studied in the paper were drawn from the literature, and
the proper citations are given in the paper. You can find much more detailed
information, including links to download many of the data sets,
here.

**A note about method tutorials**

We do not currently have any tutorial information for installing or using
these methods, beyond what we describe in the paper and what is contained in
the help files that go with the Matlab and R files themselves. That being
said, the InterSciWiki at UC Irvine has
a
good overview tutorial page that may be of some use, and Willy Lai has
created a
nice page, with R code, that works through several examples.

**Updates**

*6 December 2012*: added a link to the R package by Colin Gillespie.

*30 November 2012*: replaced plpva.r with an updated version by Neal Walfield.

*2 August 2012*: WARNING, the zeta function implementation used here is unstable for large alpha (>7) (thanks to David Gleich for pointing this out). If you need it for this range, consider using a better library function for the Hurwitz Zeta function.

*17 January 2012*: fixed a minor bug in the way plfit, plvar, plpva parse the nowarn and nosmall arguments.

*24 August 2011*: posted updated version of plfit.r, at the request of
its author Laurent Dubroca.

*4 August 2011*: posted Joel Ornstein's Python ports of plfit, plvar, plpva and plplot.

*8 October 2010*: replaced plfit.r with new version, at request of its author Laurent Dubroca.

*24 January 2010*: fixed a minor bug in how plfit.m reports the log-likelihood of the fitted data, for the discrete case, after the selection of xmin is done; posted updated version of Wim Otte's code with the same fix.

*27 November 2009*: posted Wim Otte's C++ implementation of plfit and plvar.

*1 October 2009*: added the option in plfit, plvar and plpva to 'lock' xmin to a specific value (thanks to Paul Willems for the suggestion).

*27 August 2009*: posted Adam Ginsburg's Python implementation of plfit.

*13 August 2009*: created a new page with detailed information
about obtaining copies of the 24 empirical data sets we studied.

*12 August 2009*: fixed a minor bug in the R version of plfit that would cause its
results to disagree slightly with the results from the matlab version (thanks to Naoki
Masuda for pointing it out).

*17 March 2009*: fixed a minor bug in the R version of plfit that would cause the
returned KS statistic to be incorrect when xmin=1 (thanks to Jeff Stuckman for
pointing it out).

*7 February 2009*: try-catch block in integer portion of plfit now defaults to
iterative version if the try block ever fails (thanks to Rajiv Das for the
suggestion).

*25 April 2008*: changed randht, plpva and plvar to only initialize the
pseudo-random number generator on the first time they are called.

*5 March 2008*: in the integer routines, plfit, plvar and plpva now automatically
switch to a slower but more memory efficient estimation routine when the vectorized
default routine fails (e.g., Out of Memory error when max(x) is extremely large).

*29 February 2008*: posted Laurent Dubroca's R implementation of plfit.

*17 February 2008*: posted the plplot.m function for plotting the fitted
power-law distributions against the empirical data.

*30 January 2008*: corrected typo in plpva when using hidden 'sample' option, and
reordered the commands for 'limit' and 'sample' throughout (thanks to Klaas Dellschaft
for suggestions).

*28 September 2007*: corrected typo in argument parsing for randht.m, significant
efficiency improvements to xmin estimation routine in plfit.m, plpva.m and plvar.m
(thanks to Jim Bagrow for suggestions).

*7 September 2007*: corrected interim reporting in plpva.m; changed plfit.m,
plvar.m and plpva.m to reshape input vector to column format, and to prevent using
continuous approximation in small-sample regime for discrete data.

*25 July 2007*: corrected a typo in plvar.m, typo in pareto.R, typo in
log-likelihood for discrete cut-off powerlaw and fixed small bug in a plotting
routine.

*29 June 2007*: corrected a typo in plpva.m, typo in pareto.R and updated
compilation instructions in discpowerexp.R.