Javascript Random Number With Weighted Probability

I'm trying to create a function with the following signature: function weightedRandom (target, probability) { // calculate weighted random number that is more likely to be near t

Solution 1:

TLDR: randomly choose between two easier distributions, also Inverse transform sampling

Combine two distributions

If you had a flat distribution, you could choose any value in your range equally. If you had a Gaussian distribution, you could choose a value near your Gaussian's mean. So consider randomly choosing to do one or other of these.

If you want the random value to be near the target t 80% of the time, and elsewhere the other 20%. Suppose that 'near' means within 2 standard deviations, and we'll take the variance to be v. So the range (t-2*v) to (t+2*v) needs to cover P(0.8).

Suppose we will randomly use either a flat distribution or a Gaussian distribution; the probability of a random value falling in a given range is then the sum of the two distributions, weighted by the bias of the distribution choice. If we choose the Gaussian, we will end up with a value within 2 std.dev. 95.45% of the time. If we take the Gaussian X% of the time, then the near probability Pn = P(t-2v to t+2v) = 0.9545*X + (1-X)(4v/r), where r is the full range, and (4v/r) is then the proportion of the flat distribution inside the range.

To get this Pn to 80%:

0.8 = 0.9545*X + (1-X)(4v/r). 

We have 2 unknowns, so if we also require a very near probability that the value is within 1 std.dev of the target 60% of the time, then

0.6 = 0.6827*X + (1-X)(2v/r). 

Rearranging for (2v/r):

(0.8 - 0.9545*X)/(1-X)*2 = (2v/r)
(0.6 - 0.6826*x)/(1-X) = (2v/r)

Equating and simplifying

X = 0.81546

Thus:

varrange = [0, 10];
var target = 7.0;
var stddev = 1.0;
var takeGauss = (Math.random() < 0.81546);
if(takeGauss) {
  // perform gaussian sampling (normRand has mean 0), resample if outside range
  while(1) {
    var sample = ((normRand()*stddev) + target);
    if(sample >= range[0] && sample <= range[1]) {
      return sample;
    }
  }
} else {
  // perform flat samplingreturnrange[0]+(Math.random()*(range[1]-range[0]));
}

I think this gives you the required shape, lets you choose two probabilities for near and very near probabilities, but avoids too much complexity.

As I'm being asked to provide more implementation, I've found a normal variate generator (thanks Prof Ian Neath):

functionnormRand() {
    var x1, x2, rad;

    do {
        x1 = 2 * Math.random() - 1;
        x2 = 2 * Math.random() - 1;
        rad = x1 * x1 + x2 * x2;
    } while(rad >= 1 || rad == 0);

    var c = Math.sqrt(-2 * Math.log(rad) / rad);

    return x1 * c;
};

Inverse transform sampling

The first way I considered was to use Inverse transform sampling, which I'll attempt to explain here.

Suppose we have a distribution where values from 0 to 4 are equally likely, but only half as likely as values from 4 to 10. The total probability is 4a + 6(2*a) = 1, so a=1/16:

Suppose you had a function which, when given a value between 0 and 1 in, produced a value between 0 and 10; it is still monotonic (no minima/maxima), but if you fed it every 0.01 increment from 0 to 1 you'd get a ratio of 4:6*2 = 1:3, so 3x as many values over 4 as under it. That function would look like this:

We have linear segment from z=0 to z=1/3, where x(1/3) = 4, and then a linear segment from z=1/3 to z=1 continuing on to x(1)=10. If we choose a random number z from the flat probability distribution between 0 and 1, then x(z) will be distributed with the first 1/3 of the range giving values up to 4, as required, and the remainder above.

z(x), then, is the inverse transform that takes a flat distribution and makes yields values from the desired distribution. If you want to plot it, it is x<(1/3) ? 9*x : 12*x -1.

The game is then to construct a distribution that you are happy with, and invert it to get the inverse transform, either by using pieces as above or by analytically inverting it, or some approximation (Gaussian inverse cannot be written down analytically). With this you can transform any flat-distributed sample into the desired distribution.

Sampling from the above step distribution would be done like this:

// transform 0-1 flat to 0-10 steppedfunctionstepInvTransform(z) {
    return (3*z < 1 ? 9*z : (12*z - 1));
}

// sample via inv transformvar sample = stepInvTransform(Math.random());

Solution 2:

Well. the probability is a quite easy thing to do:

//get the probability factor:
var prob = +(probability.charAt(2));//0.[2] <-- coerces 2to number

This is the max number of tries you're going to give the RNG. But then, to see how close you are getting to the target must be taken into account, which we might do with a while loop.

var targetDigit = (''+ target).charAt(2),
    prob = +((''+probability).charAt(2)),
    rn,
    closest = 0;//keep track of closest hitif (probability == 0)
    returnMath.random();
do
{
    rn = Math.random();
    //reduce chances for too great of distances, regenerate, but don't count!if (Math.abs(rn.toFixed(1).charAt(2)-targetDigit) >= target/2)
        rn = Math.random();
    if (rn.toFixed(1).charAt(2) === targetDigit)
        return rn;//on the moneyif (Math.abs(target - rn) < Math.abs(target - closest))
        closest = rn;
} while (prob--);//while probability leftreturn closest;//return the closest generated value

When I tried this code, with the parameters .8 as target, and .4 as probability, I got the following return values:

0.8257349738919855// hit!0.7696360136541552// -0.1, close0.7542727420001457// -0.1, close0.9837384401271013// +0.1, close0.3078854698178465// -0.5, far0.8178311114230021// hit!0.6079441227457084// -0.2, fairly close

Which, I take it, is close to what you'd expect.

Solution 3:

I would make an expansion of the range (Math.round()*2) and then use some other function f: [0,2] -> [0,1] that is not 1-1 mapping back to the 0-1 interval but maps more values to the area of 0.8. E.g., f(x) = c*x^2 - d*x^3 and tuning c and d such that the distribution suits you.

See this. You can just play with the values or derive the function and see how c and d interact. Of course you can use other functions for more precision like polynomials with more degrees or mix polynomials with exponential functions to get all sorts of behavior, but this is a question of what you really want, and how you would like it to behave.