The issue raised in the Geller paper is not about the ranking of the data, but rather that the Kolmogorov-Smirnov test has different critical values when the parameters of the distribution you are comparing against must be computed from data.
For an example, consider two problems: say we have some data $\mathbf{X}=(X_1,\dots,X_n)$, and we want to test (1) whether this is distributed Exp$(1)$ or (2) whether this is exponentially distributed at all. Let $F_\lambda$ denote the Exp$(\lambda)$ cdf and let $\mathbb{F}(x)=\frac1n\sum_{i=1}^n \mathbb{1}_{(-\infty,x]}(X_i)$ denote the empirical cdf. Then, to handle testing problem (1), we can use the Kolmogorov-Smirnov statistic $$
\text{KS}=\sup_x\,\lvert\mathbb{F}(x)-F_1(x)\rvert.
$$
There's no issue there.
Consider now (2). It's harder to test this; what cdf should take the role of $F_1$ when we don't know the parameter of the exponential distribution? Well, one strategy is to guess the best parameter $\lambda$, say by the MLE $\hat\lambda$ and compare against $F_{\hat\lambda}$:
$$
\text{KS}'=\sup_x\,\lvert\mathbb{F}(x)-F_{\hat\lambda}(x)\rvert.
$$
But now, this will change the distribution of the K-S statistic, and a correction must be made. Geller refers to some work by Durbin describing appropriate critical values for $\text{KS}'$ obtained from its asymptotic distribution.
Okay, so, we can see the issue for the power law paper: Clauset et al. want to estimate parameters and do a KS test. So, they have the same problem, they are not comparing against a fixed cdf but rather one that is estimated from data.
They get around this by simulating the $p$ value, rather than relying on any asymptotic results. Their procedure is as follows. Let $\mathbf{X}=(X_1,\dots,X_n)$ again be your data.
- Estimate the parameters $\hat\theta$ of the power law distribution from your data.
- Draw many simulated datasets $\mathbf{X}'_{b}$, $b=1,\dots,B$ from the power law distribution with parameters $\hat\theta$.
- For each $b=1,\dots,B$, estimate the parameters of the power law from the simulated dataset $\mathbf{X}'_{b}$, giving a new estimate $\hat\theta'_b$.
- Now, compute the KS statistic for the original dataset: $$
\text{KS}=\sup_x\,\lvert\mathbb{F}(x)-G_{\hat\theta}(x)\rvert,
$$
where $G_\theta$ is the cdf of the power law with parameters $\theta$. (Actually, in the paper it is a little more subtle, since they also model parts of the data before it enters its power law tail; I'm simplifying to get at the gist, but $G$ should be the cdf from their fitted model.)
- For each simulated dataset $b=1,\dots,B$, compute $$
\text{KS}'_b=\sup_x\,\lvert\mathbb{F}'_b(x)-G_{\hat\theta'_b}(x)\rvert,
$$
where $\mathbb{F}'_b$ is the empirical cdf of $\mathbf{X}'_b$.
- Now let $\hat{p}$ be the fraction of simulated datasets $b$ for which $\text{KS}'_b>\text{KS}$.
This is an elaboration of the procedure described on page 677 in the paper. Simulation based tests like this don't rely on any asymptotic results about the underlying statistic, so this procedure is another way of dealing with the problem (like Durbin did for the exponential -- one could also have done the same simulation procedure there, drawing simulated Exp$(\hat\lambda)$ replicates). Clauset et al. address the issue pointed out in this question and by Geller in footnote 7 at the bottom of page 677, but the point is: simulation means we don't need to worry.
I'll end with a general note: simulating $p$ values is a common way to test goodness of fit without having to do a bunch of math. The idea is, if you have some statistic, you can just see if it looks weird compared to the distribution of that statistic under your model. You just approximate that distribution by sampling from your model. It could be K-S statistic or really anything. If it looks weird, you can be like, hey my model is wrong, the empirical statistic is not what the model would like to produce. I like to plot the empirical statistic on top of the histogram of simulated ones rather than compute a $p$ value but to each their own.
