distributionspython
The Soliton distribution is a discrete probability distribution over a set $\{1,\dots, N\}$ with the probability mass function
$$
p(1)=\frac{1}{N},\qquad p(k)=\frac{1}{k(k-1)}\quad\text{for }k\in\{2,\dots, N\}
$$
I'd like to use it as part of an implementation of an LT code, ideally in Python where a uniform random number generator is available.
One of the best algorithms for sampling from a discrete distribution is the alias method.
The alias method (efficiently) precomputes a two-dimensional data structure to partition a rectangle into areas proportional to the probabilities.
In this schematic from the referenced site, a rectangle of unit height has been partitioned into four kinds of regions--as differentiated by color--in the proportions $1/2$, $1/3$, $1/12$, and $1/12$, in order to sample repeatedly from a discrete distribution with these probabilities. The vertical strips have a constant (unit) width. Each is divided into just one or two pieces. The identities of the pieces and the locations of the vertical divisions are stored in tables accessible via the column index.
The table can be sampled in two simple steps (one for each coordinate) requiring generating just two independent uniform values and $O(1)$ calculation. This improves on the $O(\log(n))$ computation needed to invert the discrete CDF as described in other replies here.
What the commenters are saying is this:
If $T$ is your translated Weibull, then you can generate variates $t_i$ from the translated Weibull with the following: $$t_i = \theta + \lambda(-\ln(1-x_i))^{1/k}\ ,$$ where the $x_i$ are your uniform random numbers.
Best Answer
If we start at $k=2$, the sums telescope, giving $1-1/k$ for the (modified) CDF. Inverting this, and taking care of the special case $k=1$, gives the following algorithm (coded in
R, I'm afraid, but you can take it as pseudocode for a Python implementation):As an example of its use (and a test), let's draw $10^5$ values for $N=10$: