# Creating PowerLaw Data

qpowerlaw[p_, alpha_, xmin_, lower_, logp_] := {
lower = True;
logp = False;
If[lower == False, p = -1 - p];
If[logp, -Exp[p]];
xmin*(1 - p)^(-1/(alpha - 1));
}

rpowerlaw[n_, xmin_, alpha_] := {
powList = Table[qpowerlaw[RandomVariate[UniformDistribution[{0,1}]]], {i,0, n}];
}


I am trying to use these two functions to properly write these lines of code in R, but I can't figure out how quite to do it.

Here is the code in R:

qpowerlaw <- function(p, alpha=2, xmin=1, lower.tail=T, log.p = F) {
if (!lower.tail)
p <- 1-p
if (log.p)
p <- exp(p)
xmin * ((1 - p) ^ (-1 / (alpha - 1)))
}
rpowerlaw <- function(n, alpha=2, xmin=1) {
qpowerlaw(runif(n, 0, 1), alpha, xmin)


}

I want to be able to test alpha = 1.8 and xmin = 10, but when I run my code. I am getting results in powlist that are less than one, which can not be right. I tried to use ParetoDistribution to properly create the data, but I could not produce the kind of results i wanted (the alpha values between the ParetoDistribution and PowerLaw are very different). I also tried PowerDistribution, which lacks the xmin parameter, so I couldn't use that function either. If anyone knows, something I could try or any errors in my code. It would be very much appreciated.

• How are you calling your function? One obvious problem I see is that qpowerlaw[] requires 5 variables, but in rpowerlaw[] you are calling it with just one. – Feyre Dec 5 '16 at 15:28
• Thank you, I didn't notice that, and it fixed up the issue. – Anthony Fasano Dec 5 '16 at 20:29

You need to check the typos in your code. (Sometimes you are translating "<-" in R to "=-" in Mathematica when it should translate just to "=".) Also, the R code results in a vector of length n rather than n+1 which is what you get from {i,0,n}.

Does the following do what you want?

qpowerlaw[p_, alpha_: 2, xmin_: 1, lower_: True, logp_: False] := {
If[! lower, p = 1 - p];
If[logp, p = Exp[p]];
xmin*(1 - p)^(-1/(alpha - 1))}

rpowerlaw[n_, alpha_: 2, xmin_: 1] :=
Table[qpowerlaw[RandomVariate[UniformDistribution[{0, 1}]], alpha, xmin][[1]], {i, n}]

rpowerlaw[12, 2, 10]
(* {322.81124530346244,48.06374170498156,18.70766795287194,18.129813319193655,
10.55503923120391,20.328561054608628,10.664461248009385,51.58041120231046,
10.103780695831814,24.197428848689704,44.91040965366457,18.41061541949269} *)