有没有一种有效的方法可以在具有给定总和或平均值的范围内生成 N 个随机整数? [英] Is there an efficient way to generate N random integers in a range that have a given sum or average?

问题描述

有没有一种有效的方法来生成 N 个整数的随机组合，使得—

Is there an efficient way to generate a random combination of N integers such that—

• 每个整数都在区间 [min, max],
• 整数的总和为 sum,
• 整数可以以任何顺序出现(例如，随机顺序)，并且
• 从满足其他要求的所有组合中随机均匀地选择该组合?

是否有类似的随机组合算法，其中整数必须按其值(而不是任何顺序)按排序顺序出现?

Is there a similar algorithm for random combinations in which the integers must appear in sorted order by their values (rather than in any order)?

(选择一个均值为mean的合适组合是一个特例，如果sum = N * mean.这个问题相当于生成一个均匀随机分区sum 分成 N 个部分，每个部分都在区间 [min, max] 中，并以任意顺序或按其值排序的顺序出现，如情况可能是这样.)

(Choosing an appropriate combination with a mean of mean is a special case, if sum = N * mean. This problem is equivalent to generating a uniform random partition of sum into N parts that are each in the interval [min, max] and appear in any order or in sorted order by their values, as the case may be.)

我知道对于以随机顺序出现的组合，可以通过以下方式解决此问题(编辑 [Apr. 27]:算法修改.):

I am aware that this problem can be solved in the following way for combinations that appear in random order (EDIT [Apr. 27]: Algorithm modified.):

1. 如果 N * max 或 N * min >sum，没有解决办法.

1. If N * max < sum or N * min > sum, there is no solution.

如果N * max == sum，则只有一种解，其中所有N个数等于max.如果N * min == sum，则只有一种解，其中所有N的数字都等于min.

If N * max == sum, there is only one solution, in which all N numbers are equal to max. If N * min == sum, there is only one solution, in which all N numbers are equal to min.

使用 Smith 和 Tromble 中给出的算法(从Unit Simplex", 2004) 生成 N 个随机非负整数，其总和为 sum - N * min.

Use the algorithm given in Smith and Tromble ("Sampling from the Unit Simplex", 2004) to generate N random non-negative integers with the sum sum - N * min.

将 min 添加到以这种方式生成的每个数字中.

Add min to each number generated this way.

如果有任何数字大于 max，请转到步骤 3.

If any number is greater than max, go to step 3.

但是，如果max 远小于sum，则该算法很慢.例如，根据我的测试(上面涉及 mean 的特殊情况的实现)，该算法平均拒绝—

However, this algorithm is slow if max is much less than sum. For example, according to my tests (with an implementation of the special case above involving mean), the algorithm rejects, on average—

• 如果 N = 7，min = 3，max = 10，sum = 42，则大约 1.6 个样本，但是
• 如果 N = 20，min = 3，max = 10，sum = 120，则大约 30.6 个样本.

• about 1.6 samples if N = 7, min = 3, max = 10, sum = 42, but
• about 30.6 samples if N = 20, min = 3, max = 10, sum = 120.

有没有办法修改这个算法，使其对大 N 有效，同时仍然满足上述要求?

Is there a way to modify this algorithm to be efficient for large N while still meeting the requirements above?

作为评论中建议的替代方法，产生有效随机组合(满足除最后一个要求之外的所有要求)的有效方法是:

As an alternative suggested in the comments, an efficient way of producing a valid random combination (that satisfies all but the last requirement) is:

1. 计算X，给定sum、min 和max 可能的有效组合数.莉>
2. 选择Y，[0, X)中的一个统一随机整数.
3. 将(unrank")Y 转换为有效的组合.

1. Calculate X, the number of valid combinations possible given sum, min, and max.
2. Choose Y, a uniform random integer in [0, X).
3. Convert ("unrank") Y to a valid combination.

但是，是否有计算有效组合(或排列)数量的公式，有没有办法将整数转换为有效组合?.

However, is there a formula for calculating the number of valid combinations (or permutations), and is there a way to convert an integer to a valid combination? .

编辑(4 月 27 日):

EDIT (Apr. 27):

阅读 Devroye 的 非均匀随机变量生成(1986 年))，我可以确认这是生成随机分区的问题.此外，第 661 页的练习 2(尤其是 E 部分)与此问题相关.

After reading Devroye's Non-Uniform Random Variate Generation (1986), I can confirm that this is a problem of generating a random partition. Also, Exercise 2 (especially part E) on page 661 is relevant to this question.

编辑(4 月 28 日):

EDIT (Apr. 28):

事实证明，我给出的算法是统一的，其中涉及的整数以随机顺序给出，而不是按值排序.由于这两个问题都具有普遍意义，因此我修改了这个问题以寻求两个问题的规范答案.

As it turned out the algorithm I gave is uniform where the integers involved are given in random order, as opposed to sorted order by their values. Since both problems are of general interest, I have modified this question to seek a canonical answer for both problems.

以下 Ruby 代码可用于验证一致性的潜在解决方案(其中 algorithm(...) 是候选算法):

The following Ruby code can be used to verify potential solutions for uniformity (where algorithm(...) is the candidate algorithm):

combos={}
permus={}
mn=0
mx=6
sum=12
for x in mn..mx
  for y in mn..mx
    for z in mn..mx
      if x+y+z==sum
        permus[[x,y,z]]=0
      end
      if x+y+z==sum and x<=y and y<=z
        combos[[x,y,z]]=0
      end
    end
  end
end

3000.times {|x|
 f=algorithm(3,sum,mn,mx)
 combos[f.sort]+=1
 permus[f]+=1
}
p combos
p permus

编辑(4 月 29 日):重新添加了当前实现的 Ruby 代码.

EDIT (Apr. 29): Re-added Ruby code of current implementation.

以下代码示例是在 Ruby 中给出的，但我的问题与编程语言无关:

The following code example is given in Ruby, but my question is independent of programming language:

def posintwithsum(n, total)
    raise if n <= 0 or total <=0
    ls = [0]
    ret = []
    while ls.length < n
      c = 1+rand(total-1)
      found = false
      for j in 1...ls.length
        if ls[j] == c
          found = true
          break
        end
      end
      if found == false;ls.push(c);end
    end
    ls.sort!
    ls.push(total)
    for i in 1...ls.length
       ret.push(ls[i] - ls[i - 1])
    end
    return ret
end

def integersWithSum(n, total)
 raise if n <= 0 or total <=0
 ret = posintwithsum(n, total + n)
 for i in 0...ret.length
    ret[i] = ret[i] - 1
 end
 return ret
end

# Generate 100 valid samples
mn=3
mx=10
sum=42
n=7
100.times {
 while true
    pp=integersWithSum(n,sum-n*mn).map{|x| x+mn }
    if !pp.find{|x| x>mx }
      p pp; break # Output the sample and break
    end
 end
}

推荐答案

这是我的 Java 解决方案.它功能齐全，包含两个生成器:用于未排序分区的 PermutationPartitionGenerator 和用于排序分区的 CombinationPartitionGenerator.您的生成器也在 SmithTromblePartitionGenerator 类中实现以进行比较.SequentialEnumerator 类按顺序枚举所有可能的分区(未排序或排序，取决于参数).我为所有这些生成器添加了全面的测试(包括您的测试用例).在大多数情况下，实现是不言自明的.如果你有任何问题，我会在几天内回答.

Here's my solution in Java. It is fully functional and contains two generators: PermutationPartitionGenerator for unsorted partitions and CombinationPartitionGenerator for sorted partitions. Your generator also implemented in the class SmithTromblePartitionGenerator for comparison. The class SequentialEnumerator enumerates all possible partitions (unsorted or sorted, depending on the parameter) in sequential order. I have added thorough tests (including your test cases) for all of these generators. The implementation is self-explainable for the most part. If you have any questions, I will answer them in couple of days.

import java.util.Random;
import java.util.function.Supplier;

public abstract class PartitionGenerator implements Supplier<int[]>{
    public static final Random rand = new Random();
    protected final int numberCount;
    protected final int min;
    protected final int range;
    protected final int sum; // shifted sum
    protected final boolean sorted;

    protected PartitionGenerator(int numberCount, int min, int max, int sum, boolean sorted) {
        if (numberCount <= 0)
            throw new IllegalArgumentException("Number count should be positive");
        this.numberCount = numberCount;
        this.min = min;
        range = max - min;
        if (range < 0)
            throw new IllegalArgumentException("min > max");
        sum -= numberCount * min;
        if (sum < 0)
            throw new IllegalArgumentException("Sum is too small");
        if (numberCount * range < sum)
            throw new IllegalArgumentException("Sum is too large");
        this.sum = sum;
        this.sorted = sorted;
    }

    // Whether this generator returns sorted arrays (i.e. combinations)
    public final boolean isSorted() {
        return sorted;
    }

    public interface GeneratorFactory {
        PartitionGenerator create(int numberCount, int min, int max, int sum);
    }
}

import java.math.BigInteger;

// Permutations with repetition (i.e. unsorted vectors) with given sum
public class PermutationPartitionGenerator extends PartitionGenerator {
    private final double[][] distributionTable;

    public PermutationPartitionGenerator(int numberCount, int min, int max, int sum) {
        super(numberCount, min, max, sum, false);
        distributionTable = calculateSolutionCountTable();
    }

    private double[][] calculateSolutionCountTable() {
        double[][] table = new double[numberCount + 1][sum + 1];
        BigInteger[] a = new BigInteger[sum + 1];
        BigInteger[] b = new BigInteger[sum + 1];
        for (int i = 1; i <= sum; i++)
            a[i] = BigInteger.ZERO;
        a[0] = BigInteger.ONE;
        table[0][0] = 1.0;
        for (int n = 1; n <= numberCount; n++) {
            double[] t = table[n];
            for (int s = 0; s <= sum; s++) {
                BigInteger z = BigInteger.ZERO;
                for (int i = Math.max(0, s - range); i <= s; i++)
                    z = z.add(a[i]);
                b[s] = z;
                t[s] = z.doubleValue();
            }
            // swap a and b
            BigInteger[] c = b;
            b = a;
            a = c;
        }
        return table;
    }

    @Override
    public int[] get() {
        int[] p = new int[numberCount];
        int s = sum; // current sum
        for (int i = numberCount - 1; i >= 0; i--) {
            double t = rand.nextDouble() * distributionTable[i + 1][s];
            double[] tableRow = distributionTable[i];
            int oldSum = s;
            // lowerBound is introduced only for safety, it shouldn't be crossed 
            int lowerBound = s - range;
            if (lowerBound < 0)
                lowerBound = 0;
            s++;
            do
                t -= tableRow[--s];
            // s can be equal to lowerBound here with t > 0 only due to imprecise subtraction
            while (t > 0 && s > lowerBound);
            p[i] = min + (oldSum - s);
        }
        assert s == 0;
        return p;
    }

    public static final GeneratorFactory factory = (numberCount, min, max,sum) ->
        new PermutationPartitionGenerator(numberCount, min, max, sum);
}

import java.math.BigInteger;

// Combinations with repetition (i.e. sorted vectors) with given sum 
public class CombinationPartitionGenerator extends PartitionGenerator {
    private final double[][][] distributionTable;

    public CombinationPartitionGenerator(int numberCount, int min, int max, int sum) {
        super(numberCount, min, max, sum, true);
        distributionTable = calculateSolutionCountTable();
    }

    private double[][][] calculateSolutionCountTable() {
        double[][][] table = new double[numberCount + 1][range + 1][sum + 1];
        BigInteger[][] a = new BigInteger[range + 1][sum + 1];
        BigInteger[][] b = new BigInteger[range + 1][sum + 1];
        double[][] t = table[0];
        for (int m = 0; m <= range; m++) {
            a[m][0] = BigInteger.ONE;
            t[m][0] = 1.0;
            for (int s = 1; s <= sum; s++) {
                a[m][s] = BigInteger.ZERO;
                t[m][s] = 0.0;
            }
        }
        for (int n = 1; n <= numberCount; n++) {
            t = table[n];
            for (int m = 0; m <= range; m++)
                for (int s = 0; s <= sum; s++) {
                    BigInteger z;
                    if (m == 0)
                        z = a[0][s];
                    else {
                        z = b[m - 1][s];
                        if (m <= s)
                            z = z.add(a[m][s - m]);
                    }
                    b[m][s] = z;
                    t[m][s] = z.doubleValue();
                }
            // swap a and b
            BigInteger[][] c = b;
            b = a;
            a = c;
        }
        return table;
    }

    @Override
    public int[] get() {
        int[] p = new int[numberCount];
        int m = range; // current max
        int s = sum; // current sum
        for (int i = numberCount - 1; i >= 0; i--) {
            double t = rand.nextDouble() * distributionTable[i + 1][m][s];
            double[][] tableCut = distributionTable[i];
            if (s < m)
                m = s;
            s -= m;
            while (true) {
                t -= tableCut[m][s];
                // m can be 0 here with t > 0 only due to imprecise subtraction
                if (t <= 0 || m == 0)
                    break;
                m--;
                s++;
            }
            p[i] = min + m;
        }
        assert s == 0;
        return p;
    }

    public static final GeneratorFactory factory = (numberCount, min, max, sum) ->
        new CombinationPartitionGenerator(numberCount, min, max, sum);
}

import java.util.*;

public class SmithTromblePartitionGenerator extends PartitionGenerator {
    public SmithTromblePartitionGenerator(int numberCount, int min, int max, int sum) {
        super(numberCount, min, max, sum, false);
    }

    @Override
    public int[] get() {
        List<Integer> ls = new ArrayList<>(numberCount + 1);
        int[] ret = new int[numberCount];
        int increasedSum = sum + numberCount;
        while (true) {
            ls.add(0);
            while (ls.size() < numberCount) {
                int c = 1 + rand.nextInt(increasedSum - 1);
                if (!ls.contains(c))
                    ls.add(c);
            }
            Collections.sort(ls);
            ls.add(increasedSum);
            boolean good = true;
            for (int i = 0; i < numberCount; i++) {
                int x = ls.get(i + 1) - ls.get(i) - 1;
                if (x > range) {
                    good = false;
                    break;
                }
                ret[i] = x;
            }
            if (good) {
                for (int i = 0; i < numberCount; i++)
                    ret[i] += min;
                return ret;
            }
            ls.clear();
        }
    }

    public static final GeneratorFactory factory = (numberCount, min, max, sum) ->
        new SmithTromblePartitionGenerator(numberCount, min, max, sum);
}

import java.util.Arrays;

// Enumerates all partitions with given parameters
public class SequentialEnumerator extends PartitionGenerator {
    private final int max;
    private final int[] p;
    private boolean finished;

    public SequentialEnumerator(int numberCount, int min, int max, int sum, boolean sorted) {
        super(numberCount, min, max, sum, sorted);
        this.max = max;
        p = new int[numberCount];
        startOver();
    }

    private void startOver() {
        finished = false;
        int unshiftedSum = sum + numberCount * min;
        fillMinimal(0, Math.max(min, unshiftedSum - (numberCount - 1) * max), unshiftedSum);
    }

    private void fillMinimal(int beginIndex, int minValue, int fillSum) {
        int fillRange = max - minValue;
        if (fillRange == 0)
            Arrays.fill(p, beginIndex, numberCount, max);
        else {
            int fillCount = numberCount - beginIndex;
            fillSum -= fillCount * minValue;
            int maxCount = fillSum / fillRange;
            int maxStartIndex = numberCount - maxCount;
            Arrays.fill(p, maxStartIndex, numberCount, max);
            fillSum -= maxCount * fillRange;
            Arrays.fill(p, beginIndex, maxStartIndex, minValue);
            if (fillSum != 0)
                p[maxStartIndex - 1] = minValue + fillSum;
        }
    }

    @Override
    public int[] get() { // returns null when there is no more partition, then starts over
        if (finished) {
            startOver();
            return null;
        }
        int[] pCopy = p.clone();
        if (numberCount > 1) {
            int i = numberCount;
            int s = p[--i];
            while (i > 0) {
                int x = p[--i];
                if (x == max) {
                    s += x;
                    continue;
                }
                x++;
                s--;
                int minRest = sorted ? x : min;
                if (s < minRest * (numberCount - i - 1)) {
                    s += x;
                    continue;
                }
                p[i++]++;
                fillMinimal(i, minRest, s);
                return pCopy;
            }
        }
        finished = true;
        return pCopy;
    }

    public static final GeneratorFactory permutationFactory = (numberCount, min, max, sum) ->
        new SequentialEnumerator(numberCount, min, max, sum, false);
    public static final GeneratorFactory combinationFactory = (numberCount, min, max, sum) ->
        new SequentialEnumerator(numberCount, min, max, sum, true);
}

import java.util.*;
import java.util.function.BiConsumer;
import PartitionGenerator.GeneratorFactory;

public class Test {
    private final int numberCount;
    private final int min;
    private final int max;
    private final int sum;
    private final int repeatCount;
    private final BiConsumer<PartitionGenerator, Test> procedure;

    public Test(int numberCount, int min, int max, int sum, int repeatCount,
            BiConsumer<PartitionGenerator, Test> procedure) {
        this.numberCount = numberCount;
        this.min = min;
        this.max = max;
        this.sum = sum;
        this.repeatCount = repeatCount;
        this.procedure = procedure;
    }

    @Override
    public String toString() {
        return String.format("=== %d numbers from [%d, %d] with sum %d, %d iterations ===",
                numberCount, min, max, sum, repeatCount);
    }

    private static class GeneratedVector {
        final int[] v;

        GeneratedVector(int[] vect) {
            v = vect;
        }

        @Override
        public int hashCode() {
            return Arrays.hashCode(v);
        }

        @Override
        public boolean equals(Object obj) {
            if (this == obj)
                return true;
            return Arrays.equals(v, ((GeneratedVector)obj).v);
        }

        @Override
        public String toString() {
            return Arrays.toString(v);
        }
    }

    private static final Comparator<Map.Entry<GeneratedVector, Integer>> lexicographical = (e1, e2) -> {
        int[] v1 = e1.getKey().v;
        int[] v2 = e2.getKey().v;
        int len = v1.length;
        int d = len - v2.length;
        if (d != 0)
            return d;
        for (int i = 0; i < len; i++) {
            d = v1[i] - v2[i];
            if (d != 0)
                return d;
        }
        return 0;
    };

    private static final Comparator<Map.Entry<GeneratedVector, Integer>> byCount =
            Comparator.<Map.Entry<GeneratedVector, Integer>>comparingInt(Map.Entry::getValue)
            .thenComparing(lexicographical);

    public static int SHOW_MISSING_LIMIT = 10;

    private static void checkMissingPartitions(Map<GeneratedVector, Integer> map, PartitionGenerator reference) {
        int missingCount = 0;
        while (true) {
            int[] v = reference.get();
            if (v == null)
                break;
            GeneratedVector gv = new GeneratedVector(v);
            if (!map.containsKey(gv)) {
                if (missingCount == 0)
                    System.out.println(" Missing:");
                if (++missingCount > SHOW_MISSING_LIMIT) {
                    System.out.println("  . . .");
                    break;
                }
                System.out.println(gv);
            }
        }
    }

    public static final BiConsumer<PartitionGenerator, Test> distributionTest(boolean sortByCount) {
        return (PartitionGenerator gen, Test test) -> {
            System.out.print("
" + getName(gen) + "

");
            Map<GeneratedVector, Integer> combos = new HashMap<>();
            // There's no point of checking permus for sorted generators
            // because they are the same as combos for them
            Map<GeneratedVector, Integer> permus = gen.isSorted() ? null : new HashMap<>();
            for (int i = 0; i < test.repeatCount; i++) {
                int[] v = gen.get();
                if (v == null && gen instanceof SequentialEnumerator)
                    break;
                if (permus != null) {
                    permus.merge(new GeneratedVector(v), 1, Integer::sum);
                    v = v.clone();
                    Arrays.sort(v);
                }
                combos.merge(new GeneratedVector(v), 1, Integer::sum);
            }
            Set<Map.Entry<GeneratedVector, Integer>> sortedEntries = new TreeSet<>(
                    sortByCount ? byCount : lexicographical);
            System.out.println("Combos" + (gen.isSorted() ? ":" : " (don't have to be uniform):"));
            sortedEntries.addAll(combos.entrySet());
            for (Map.Entry<GeneratedVector, Integer> e : sortedEntries)
                System.out.println(e);
            checkMissingPartitions(combos, test.getGenerator(SequentialEnumerator.combinationFactory));
            if (permus != null) {
                System.out.println("
Permus:");
                sortedEntries.clear();
                sortedEntries.addAll(permus.entrySet());
                for (Map.Entry<GeneratedVector, Integer> e : sortedEntries)
                    System.out.println(e);
                checkMissingPartitions(permus, test.getGenerator(SequentialEnumerator.permutationFactory));
            }
        };
    }

    public static final BiConsumer<PartitionGenerator, Test> correctnessTest =
        (PartitionGenerator gen, Test test) -> {
        String genName = getName(gen);
        for (int i = 0; i < test.repeatCount; i++) {
            int[] v = gen.get();
            if (v == null && gen instanceof SequentialEnumerator)
                v = gen.get();
            if (v.length != test.numberCount)
                throw new RuntimeException(genName + ": array of wrong length");
            int s = 0;
            if (gen.isSorted()) {
                if (v[0] < test.min || v[v.length - 1] > test.max)
                    throw new RuntimeException(genName + ": generated number is out of range");
                int prev = test.min;
                for (int x : v) {
                    if (x < prev)
                        throw new RuntimeException(genName + ": unsorted array");
                    s += x;
                    prev = x;
                }
            } else
                for (int x : v) {
                    if (x < test.min || x > test.max)
                        throw new RuntimeException(genName + ": generated number is out of range");
                    s += x;
                }
            if (s != test.sum)
                throw new RuntimeException(genName + ": wrong sum");
        }
        System.out.format("%30s :   correctness test passed%n", genName);
    };

    public static final BiConsumer<PartitionGenerator, Test> performanceTest =
        (PartitionGenerator gen, Test test) -> {
        long time = System.nanoTime();
        for (int i = 0; i < test.repeatCount; i++)
            gen.get();
        time = System.nanoTime() - time;
        System.out.format("%30s : %8.3f s %10.0f ns/test%n", getName(gen), time * 1e-9, time * 1.0 / test.repeatCount);
    };

    public PartitionGenerator getGenerator(GeneratorFactory factory) {
        return factory.create(numberCount, min, max, sum);
    }

    public static String getName(PartitionGenerator gen) {
        String name = gen.getClass().getSimpleName();
        if (gen instanceof SequentialEnumerator)
            return (gen.isSorted() ? "Sorted " : "Unsorted ") + name;
        else
            return name;
    }

    public static GeneratorFactory[] factories = { SmithTromblePartitionGenerator.factory,
            PermutationPartitionGenerator.factory, CombinationPartitionGenerator.factory,
            SequentialEnumerator.permutationFactory, SequentialEnumerator.combinationFactory };

    public static void main(String[] args) {
        Test[] tests = {
                            new Test(3, 0, 3, 5, 3_000, distributionTest(false)),
                            new Test(3, 0, 6, 12, 3_000, distributionTest(true)),
                            new Test(50, -10, 20, 70, 2_000, correctnessTest),
                            new Test(7, 3, 10, 42, 1_000_000, performanceTest),
                            new Test(20, 3, 10, 120, 100_000, performanceTest)
                       };
        for (Test t : tests) {
            System.out.println(t);
            for (GeneratorFactory factory : factories) {
                PartitionGenerator candidate = t.getGenerator(factory);
                t.procedure.accept(candidate, t);
            }
            System.out.println();
        }
    }
}

您可以在 Ideone 上试试这个.

这篇关于有没有一种有效的方法可以在具有给定总和或平均值的范围内生成 N 个随机整数?的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！