<p>There is a package named powerRlaw in R by Colin Gillespie. This package is well documented and contains a lot of example to use each function. Very straightforward.</p>

<p><a href="http://cran.r-project.org/web/packages/poweRlaw/" rel="nofollow noreferrer">http://cran.r-project.org/web/packages/poweRlaw/</a></p>

<p>For example in R as the documentation said, the following code get data from the file <strong>full_path_of_file_name</strong> and estimate xmin and alpha and get p-value as proposed by <a href="https://arxiv.org/pdf/0706.1062.pdf" rel="nofollow noreferrer">Clauset and al. (2009)</a></p>

<pre><code>library("poweRLaw")

words = read.table(&lt;full_path_of_file_name&gt;)
m_plwords = displ$new(words$V1)         # discrete power law fitting
est_plwords = estimate_xmin(m_plwords)  # get xmin and alpha

# here we have the goodness-of-fit test p-value
# as proposed by Clauset and al. (2009)
bs_p = bootstrap_p(m_plwords)              
</code></pre>
Goodness of fit test for power law distribution in R
<p>I have a network for which I fit into a power-law using igraph software:</p>

<pre><code>plf = power.law.fit(degree_dist, impelementation = "plfit")
</code></pre>

<p>The plf variable now holds the following variables:</p>

<pre><code>$continuous
[1] TRUE
$alpha
[1] 1.63975
$xmin
[1] 0.03
$logLik
[1] 4.037563
$KS.stat
[1] 0.1721117
$KS.p
[1] 0.9984284
</code></pre>

<p>The igraph manual explains these variables:</p>

<pre><code>xmin = the lower bound for fitting the power-law
alpha =  the exponent of the fitted power-law distribution
logLik =  the log-likelihood of the fitted parameters
KS.stat =  the test statistic of a Kolmogorov-Smirnov test that compares the fitted  distribution with the input vector. Smaller scores denote better fit
KS.p = the p-value of the Kolmogorov-Smirnov test. Small p-values (less than 0.05) indicate that the test rejected the hypothesis that the original data could have been drawn from the fitted power-law distribution
</code></pre>

<p>I would like to do a "goodness of fit" test on this power law fit. But I am not sure how to do this, and although I found this question already asked in online forums, it usually remains unanswered.</p>

<p>I think one way to do this would be to do a chisq.test(x,y). One input parameter (say x) would be the degree_dist variable (the observed degree distribution of the network). The other input parameter (say y) would be the fitted power law equation, which is supposed to be of form P(k) = mk^a.</p>

<p>I am not sure whether this is a sound approach, and if so, I need advice on how to construct the fitted power law equation.</p>

<p>In case it helps, the degree_dist of my network was:</p>

<pre><code> 0.00 0.73 0.11 0.05 0.02 0.02 0.03 0.02 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00        0.01 0.00 0.00 0.00 0.01
</code></pre>

<p>(These are frequencies that degrees of 0-21 occurred in the network. (For example, 73% of nodes has degree 1, 1% of nodes had degree 21).</p>

<p><em><strong></em>**<em>*</em>**<em>*</em>***<em></strong>  EDIT <strong></em>**<em>*</em>**<em>*</em>**<em>*</em>****</strong> </p>

<p>I am unsure whether it was a mistake above to use degree_dist to calculate plf. In case it is, I also ran the same function using the degrees from the 100 nodes in my network:</p>

<pre><code>plf = power.law.fit(pure_deg, impelementation = "plfit")
</code></pre>

<p>where, pure_deg is:</p>

<pre><code>  21  7  5  6 17  3  6  6  2  5  4  3  7  4  3  2  2  2  2  3  2  3  2  2  2  2  2  1  1  1  1  1  1 1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1 1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1 1
</code></pre>

<p>This leads to output of:</p>

<pre><code>$continuous
[1] FALSE
$alpha
[1] 2.362445
$xmin
[1] 1
$logLik
[1] -114.6303
$KS.stat
[1] 0.02293443
$KS.p
[1] 1
</code></pre>