¿Hay una manera eficiente de generar N enteros aleatorios en un rango que tenga una suma o promedio dado?

14

¿Hay una manera eficiente de generar una combinación aleatoria de N enteros de tal manera que:

  • cada entero está en el intervalo [ min, max],
  • los enteros tienen una suma de sum,
  • los enteros pueden aparecer en cualquier orden (p. ej., orden aleatorio) y
  • la combinación se elige de manera uniforme al azar entre todas las combinaciones que cumplen con los otros requisitos?

¿Existe un algoritmo similar para combinaciones aleatorias en el que los enteros deben aparecer en orden ordenado por sus valores (en lugar de en cualquier orden)?

(Elegir una combinación apropiada con una media de meanes un caso especial, si sum = N * mean. Este problema es equivalente a generar una partición aleatoria uniforme de sumN partes que están en el intervalo [ min, max] y aparecen en cualquier orden o en orden ordenado por su valores, según sea el caso.)

Soy consciente de que este problema se puede resolver de la siguiente manera para combinaciones que aparecen en orden aleatorio (EDIT [27 de abril]: Algoritmo modificado):

  1. Si N * max < sumo N * min > sum, no hay solución.

  2. Si N * max == sum, solo hay una solución, en la que todos los Nnúmeros son iguales max. Si N * min == sum, solo hay una solución, en la que todos los Nnúmeros son iguales min.

  3. Utilice el algoritmo proporcionado en Smith y Tromble ("Muestreo de la unidad simple", 2004) para generar N enteros no negativos aleatorios con la suma sum - N * min.

  4. Agregue mina cada número generado de esta manera.

  5. Si cualquier número es mayor que max, vaya al paso 3.

Sin embargo, este algoritmo es lento si maxes mucho menor que sum. Por ejemplo, de acuerdo con mis pruebas (con una implementación del caso especial mencionado anteriormente mean), el algoritmo rechaza, en promedio:

  • alrededor de 1.6 muestras si N = 7, min = 3, max = 10, sum = 42, pero
  • aproximadamente 30,6 muestras si N = 20, min = 3, max = 10, sum = 120.

¿Hay alguna forma de modificar este algoritmo para que sea eficiente para N grande mientras se cumplen los requisitos anteriores?

EDITAR:

Como alternativa sugerida en los comentarios, una forma eficiente de producir una combinación aleatoria válida (que satisfaga todos los requisitos excepto el último) es:

  1. Calcular X, el número de combinaciones válidas posible dada sum, miny max.
  2. Elija Y, un entero aleatorio uniforme en [0, X).
  3. Convierta ("sin carga") Yen una combinación válida.

Sin embargo, ¿existe una fórmula para calcular el número de combinaciones válidas (o permutaciones), y hay alguna forma de convertir un número entero en una combinación válida? [EDITAR (28 de abril): lo mismo para las permutaciones que para las combinaciones].

EDITAR (27 de abril):

Después de leer la generación de varianza aleatoria no uniforme de Devroye (1986), puedo confirmar que este es un problema de generar una partición aleatoria. Además, el ejercicio 2 (especialmente la parte E) en la página 661 es relevante para esta pregunta.

EDITAR (28 de abril):

Al final resultó que el algoritmo que di es uniforme, donde los enteros involucrados se dan en orden aleatorio , en lugar de ordenarlos por sus valores . Como ambos problemas son de interés general, he modificado esta pregunta para buscar una respuesta canónica para ambos problemas.

El siguiente código Ruby se puede usar para verificar posibles soluciones para la uniformidad (donde algorithm(...)está el algoritmo candidato):

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

EDITAR (29 de abril): se volvió a agregar el código Ruby de la implementación actual.

El siguiente ejemplo de código se da en Ruby, pero mi pregunta es independiente del lenguaje de programación:

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
}
Peter O.
fuente
¿Podría aclarar su tercer requisito? ¿Necesita una uniformidad entre todas las combinaciones posibles (incluidas las que tienen la media incorrecta) o entre todas las combinaciones válidas (es decir, las que tienen la media correcta)?
user58697
Todas las combinaciones válidas, es decir, todas las combinaciones que cumplen los otros requisitos.
Peter O.
Si tuviéramos una manera de contar y eliminar particiones de una suma restringida a N enteros en [min, max], ¿elegir una de esas particiones al azar y sin clasificar representaría una distribución uniforme, y sería eso más eficiente que su método actual? ¿Qué tan grande puede ser la suma y N?
עדלעד ברקן
No sé a qué te refieres con "desmarcar particiones de una suma", y no estoy al tanto de una prueba de que hacerlo resulte en una distribución uniforme en el sentido de esta pregunta. Para esta pregunta, tanto sumy Nson prácticamente ilimitada (dentro de lo razonable). Estoy buscando una respuesta canónica porque el problema subyacente aparece en muchas preguntas formuladas en Stack Overflow, incluidas esta y esta . @ גלעדברקן
Peter O.
Si le damos a cada combinación posible un "rango" (o índice) en una disposición ordenada de todos ellos, "sin clasificar" significaría generar la combinación, dado su rango (y N, min y max, por supuesto). ¿Por qué una elección de una de todas las combinaciones posibles no se ajusta a una distribución uniforme?
עדלעד ברקן

Respuestas:

5

Aquí está mi solución en Java. Es completamente funcional y contiene dos generadores: PermutationPartitionGeneratorpara particiones sin clasificar y CombinationPartitionGeneratorpara particiones ordenadas. Su generador también se implementó en la clase SmithTromblePartitionGeneratorpara comparación. La clase SequentialEnumeratorenumera todas las particiones posibles (sin clasificar u ordenadas, según el parámetro) en orden secuencial. He agregado pruebas exhaustivas (incluidos sus casos de prueba) para todos estos generadores. La implementación es autoexplicable en su mayor parte. Si tiene alguna pregunta, las responderé en un par de días.

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("\n" + getName(gen) + "\n\n");
            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("\nPermus:");
                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();
        }
    }
}

Puedes probar esto en Ideone .

John McClane
fuente
Gracias por tu respuesta; funciona bien. He descrito el generador de permutación en otra respuesta aquí; respondió otra pregunta con su ayuda; y pronto incluirá su algoritmo en el código de muestra de Python para mi artículo sobre métodos de generación aleatoria.
Peter O.
Solo para aclarar. ¿Este algoritmo se basa en generar todas las particiones / composiciones posibles para muestrear?
Joseph Wood
@JosephWood No, depende de contarlos a todos. Esto se hace solo una vez en la inicialización del generador y es bastante efectivo porque utiliza el enfoque de programación dinámica.
John McClane hace
¿Cómo puede la programación dinámica resolver el problema relacionado de elegir una partición aleatoria uniforme de 'suma' en N enteros elegidos al azar con reemplazo de una lista ( ejemplo ) o sin reemplazo ( ejemplo ), o cómo puede resolverse ese problema?
Peter O.
@PeterO. Debe contar todas las particiones posibles a través del mismo método que en mi algoritmo, pero esta vez necesita restar solo los números permitidos de la suma. Esto es demasiado largo para comentar, puede hacer una pregunta por separado. Sospecho que uno puede resolver cuatro problemas diferentes a través del mismo enfoque. Suponga que tiene una lista de enteros distintos para elegir (esto es solo un rango continuo en esta pregunta). Luego, puede generar matrices aleatorias de una longitud dada que consisten en números de esta lista con la suma dada si las matrices deben clasificarse / no clasificarse y permitir / rechazar una repetición.
John McClane hace
3

Aquí está el algoritmo de PermutationPartitionGenerator de John McClane, en otra respuesta en esta página. Tiene dos fases, a saber, una fase de configuración y una fase de muestreo, y genera nnúmeros aleatorios en [ min, max] con la suma sum, donde los números se enumeran en orden aleatorio.

Fase de configuración: Primero, se crea una tabla de solución utilizando las siguientes fórmulas ( t(y, x)donde yestá en [0, n] y xestá en [0, sum - n * min]):

  • t (0, j) = 1 si j == 0; 0 de lo contrario
  • t (i, j) = t (i-1, j) + t (i-1, j-1) + ... + t (i-1, j- (max-min))

Aquí, t (y, x) almacena la probabilidad relativa de que la suma de ynúmeros (en el rango apropiado) sea igual x. Esta probabilidad es relativa a todos los t (y, x) con el mismo y.

Fase de muestreo: aquí generamos una muestra de nnúmeros. Ajuste sa sum - n * min, luego para cada posición i, comenzando n - 1y trabajando hacia atrás a 0:

  • Establezca vun entero aleatorio en [0, t (i + 1, s)).
  • Establecer ra min.
  • Resta t (i, s) de v.
  • Mientras vpermanezca 0 o mayor, reste t (i, s-1) de v, sume 1 ry reste 1 de s.
  • El número en la posición ide la muestra se establece en r.

EDITAR:

Parece que con cambios triviales en el algoritmo anterior, es posible que cada número aleatorio use un rango separado en lugar de usar el mismo rango para todos ellos:

Cada número aleatorio en las posiciones i∈ [0, n) tiene un valor mínimo min (i) y un valor máximo max (i).

Let adjsum= sum- Σmin (i).

Fase de configuración: Primero, se crea una tabla de solución utilizando las siguientes fórmulas ( t(y, x)donde yestá en [0, n] y xestá en [0, adjsum]):

  • t (0, j) = 1 si j == 0; 0 de lo contrario
  • t (i, j) = t (i-1, j) + t (i-1, j-1) + ... + t (i-1, j- (max (i-1) -min (i -1)) )

La fase de muestreo es exactamente la misma que antes, excepto que establecemos sen adjsum(en lugar de sum - n * min) y establecemos ren min (i) (en lugar de min).


EDITAR:

Para CombinationPartitionGenerator de John McClane, las fases de configuración y muestreo son las siguientes.

Fase de configuración: Primero, se construye una tabla de solución usando las siguientes fórmulas ( t(z, y, x)donde zestá en [0, n], yestá en [0, max - min] y xestá en [0, sum - n * min]):

  • t (0, j, k) = 1 si k == 0; 0 de lo contrario
  • t (i, 0, k) = t (i - 1, 0, k)
  • t (i, j, k) = t (i, j-1, k) + t (i - 1, j, k - j)

Fase de muestreo: aquí generamos una muestra de nnúmeros. Establezca sto sum - n * miny mrangeto max - min, luego para cada posición i, comenzando n - 1y trabajando hacia atrás a 0:

  • Establecer ven un entero aleatorio en [0, t (i + 1, rango, s)).
  • Establecer mrangeen min ( mrange, s)
  • Restar mrangede s.
  • Establecer ra min + mrange.
  • T Reste ( i, mrange, s) a partir de v.
  • Mientras vrestos 0 o mayor, añadir 1 a s, restar 1 a ry 1 de mrange, a continuación, t restar ( i, mrange, s) a partir de v.
  • El número en la posición ide la muestra se establece en r.
Peter O.
fuente
2

No he probado esto, por lo que no es realmente una respuesta, solo algo para probar que es demasiado largo para caber en un comentario. Comience con una matriz que cumpla los dos primeros criterios y juegue con ella para que cumpla con los dos primeros, pero es mucho más aleatorio.

Si la media es un número entero, entonces su matriz inicial puede ser [4, 4, 4, ... 4] o tal vez [3, 4, 5, 3, 4, 5, ... 5, 8, 0] o Algo simple como eso. Para una media de 4.5, intente [4, 5, 4, 5, ... 4, 5].

A continuación, elija un par de números num1y num2, en la matriz. Probablemente, el primer número debe tomarse en orden, ya que con el shuffle de Fisher-Yates, el segundo número debe elegirse al azar. Tomar el primer número en orden asegura que cada número se elija al menos una vez.

Ahora calcule max-num1y num2-min. Esas son las distancias de los dos números a los límites maxy min. Establecer limiten la menor de las dos distancias. Ese es el cambio máximo permitido que no pondrá uno u otro número fuera de los límites permitidos. Si limites cero, omita este par.

Elija un entero aleatorio en el rango [1, limit]: llámelo change. Omito 0 del rango seleccionable ya que no tiene ningún efecto. Las pruebas pueden mostrar que obtienes una mejor aleatoriedad al incluirla; No estoy seguro.

Ahora listo num1 <- num1 + changey num2 <- num2 - change. Eso no afectará el valor medio y todos los elementos de la matriz todavía están dentro de los límites requeridos.

Deberá ejecutar la matriz completa al menos una vez. Las pruebas deben mostrar si necesita ejecutarlo más de una vez para obtener algo lo suficientemente aleatorio.

ETA: incluye pseudocódigo

// Set up the array.
resultAry <- new array size N
for (i <- 0 to N-1)
  // More complex initial setup schemes are possible here.
  resultAry[i] <- mean
rof

// Munge the array entries.
for (ix1 <- 0 to N-1)  // ix1 steps through the array in order.

  // Pick second entry different from first.
  repeat
    ix2 <- random(0, N-1)
  until (ix2 != ix1)

  // Calculate size of allowed change.
  hiLimit <- max - resultAry[ix1]
  loLimit <- resultAry[ix2] - min
  limit <- minimum(hiLimit, loLimit)
  if (limit == 0)
    // No change possible so skip.
    continue loop with next ix1
  fi

  // Change the two entries keeping same mean.
  change <- random(1, limit)  // Or (0, limit) possibly.
  resultAry[ix1] <- resultAry[ix1] + change
  resultAry[ix2] <- resultAry[ix2] - change

rof

// Check array has been sufficiently munged.
if (resultAry not random enough)
  munge the array again
fi
rossum
fuente
Lo he probado y, desafortunadamente, su algoritmo no forma una distribución uniforme de todas las soluciones, sin importar cuántas iteraciones haga.
Peter O.
Oh bien. Valió la pena intentarlo de todos modos. :(
rossum
2

Como señala el OP, la capacidad de descargar eficientemente es muy poderosa. Si podemos hacerlo, se puede generar una distribución uniforme de particiones en tres pasos (reafirmando lo que el OP ha establecido en la pregunta):

  1. Calcule el número total, M , de particiones de longitud N del número de summodo que las partes estén en el rango [ min, max].
  2. Genere una distribución uniforme de enteros a partir de [1, M].
  3. Descargue cada número entero del paso 2 en su partición respectiva.

A continuación, sólo se centran en la generación de la n º partición ya que hay una copiosa cantidad de información sobre la generación de una distribución uniforme de número entero en un rango determinado. Aquí hay un C++algoritmo de clasificación simple que debería ser fácil de traducir a otros idiomas (Nota: todavía no he descubierto cómo eliminar el caso de composición (es decir, el orden importa)).

std::vector<int> unRank(int n, int m, int myMax, int nth) {

    std::vector<int> z(m, 0);
    int count = 0;
    int j = 0;

    for (int i = 0; i < z.size(); ++i) {
        int temp = pCount(n - 1, m - 1, myMax);

        for (int r = n - m, k = myMax - 1;
             (count + temp) < nth && r > 0 && k; r -= m, --k) {

            count += temp;
            n = r;
            myMax = k;
            ++j;
            temp = pCount(n - 1, m - 1, myMax);
        }

        --m;
        --n;
        z[i] = j;
    }

    return z;
}

La pCountfunción de caballo de batalla viene dada por:

int pCount(int n, int m, int myMax) {

    if (myMax * m < n) return 0;
    if (myMax * m == n) return 1;

    if (m < 2) return m;
    if (n < m) return 0;
    if (n <= m + 1) return 1;

    int niter = n / m;
    int count = 0;

    for (; niter--; n -= m, --myMax) {
        count += pCount(n - 1, m - 1, myMax);
    }

    return count;
}

Esta función se basa en la excelente respuesta a ¿Existe un algoritmo eficiente para la partición de enteros con un número restringido de partes? por el usuario @ m69_snarky_and_unwelcoming. El que se da arriba es una ligera modificación del algoritmo simple (el que no tiene memoria). Esto se puede modificar fácilmente para incorporar la memorización para una mayor eficiencia. Dejaremos esto por ahora y nos centraremos en la parte de clasificación.

Explicación de unRank

Primero notamos que hay un mapeo uno a uno desde las particiones de longitud N del número de summanera que las partes están en el rango [ min, max] a las particiones restringidas de longitud N del número sum - N * (min - 1)con partes en [ 1, max - (min - 1)].

Como un pequeño ejemplo, considere las particiones 50de longitud 4tal que el min = 10y el max = 15. Esto tendrá la misma estructura que las particiones restringidas 50 - 4 * (10 - 1) = 14de longitud 4con la parte máxima igual a 15 - (10 - 1) = 6.

10   10   15   15   --->>    1    1    6    6
10   11   14   15   --->>    1    2    5    6
10   12   13   15   --->>    1    3    4    6
10   12   14   14   --->>    1    3    5    5
10   13   13   14   --->>    1    4    4    5
11   11   13   15   --->>    2    2    4    6
11   11   14   14   --->>    2    2    5    5
11   12   12   15   --->>    2    3    3    6
11   12   13   14   --->>    2    3    4    5
11   13   13   13   --->>    2    4    4    4
12   12   12   14   --->>    3    3    3    5
12   12   13   13   --->>    3    3    4    4

Con esto en mente, para contar fácilmente, podríamos agregar un paso 1a para traducir el problema al caso de "unidad" si lo desea.

Ahora, simplemente tenemos un problema de conteo. Como se muestra brillantemente en @ m69, el recuento de particiones se puede lograr fácilmente dividiendo el problema en problemas más pequeños. La función que proporciona @ m69 nos proporciona el 90% del camino, solo tenemos que averiguar qué hacer con la restricción adicional de que hay un límite. Aquí es donde obtenemos:

int pCount(int n, int m, int myMax) {

    if (myMax * m < n) return 0;
    if (myMax * m == n) return 1;

También tenemos que tener en cuenta que myMaxdisminuirá a medida que avancemos. Esto tiene sentido si nos fijamos en la 6 ª partición anterior:

2   2   4   6

Para contar el número de particiones de aquí en adelante, debemos seguir aplicando la traducción al caso de "unidad". Esto se ve así:

1   1   3   5

Donde como el paso anterior, teníamos un máximo de 6, ahora solo consideramos un máximo de 5.

Con esto en mente, desentrañar la partición no es diferente a descalificar una permutación o combinación estándar. Debemos poder contar el número de particiones en una sección determinada. Por ejemplo, para contar el número de particiones que comienzan con lo 10anterior, todo lo que hacemos es eliminar el 10en la primera columna:

10   10   15   15
10   11   14   15
10   12   13   15
10   12   14   14
10   13   13   14

10   15   15
11   14   15
12   13   15
12   14   14
13   13   14

Traducir a la caja de la unidad:

1   6   6
2   5   6
3   4   6
3   5   5
4   4   5

y llama pCount:

pCount(13, 3, 6) = 5

Dado un entero aleatorio para descargar, continuamos calculando el número de particiones en secciones cada vez más pequeñas (como hicimos anteriormente) hasta que hayamos llenado nuestro vector índice.

Ejemplos

Teniendo en cuenta min = 3, max = 10, n = 7, y sum = 42, aquí es una Ideone demo que genera 20 particiones aleatorias. La salida está abajo:

42: 3 3 6 7 7 8 8 
123: 4 4 6 6 6 7 9 
2: 3 3 3 4 9 10 10 
125: 4 4 6 6 7 7 8 
104: 4 4 4 6 6 8 10 
74: 3 4 6 7 7 7 8 
47: 3 4 4 5 6 10 10 
146: 5 5 5 5 6 7 9 
70: 3 4 6 6 6 7 10 
134: 4 5 5 6 6 7 9 
136: 4 5 5 6 7 7 8 
81: 3 5 5 5 8 8 8 
122: 4 4 6 6 6 6 10 
112: 4 4 5 5 6 8 10 
147: 5 5 5 5 6 8 8 
142: 4 6 6 6 6 7 7 
37: 3 3 6 6 6 9 9 
67: 3 4 5 6 8 8 8 
45: 3 4 4 4 8 9 10 
44: 3 4 4 4 7 10 10

El índice lexicográfico está a la izquierda y la partición sin clasificar a la derecha.

Joseph Wood
fuente
1
Como resultado, esta es una muy buena alternativa, y de hecho se vuelve eficiente con la memorización.
Peter O.
1
Gran observación sobre el mapeo uno a uno.
גלעד ברקן
0

Si genera 0≤a≤1 de los valores aleatorios en el rango [l, x-1] de manera uniforme, y 1-a de los valores aleatorios en el rango [x, h] de manera uniforme, la media esperada sería:

m = ((l+x-1)/2)*a + ((x+h)/2)*(1-a)

Entonces, si quieres una m específica, puedes jugar con a y x.

Por ejemplo, si establece x = m: a = (hm) / (h-l + 1).

Para garantizar una probabilidad más cercana a la uniforme para diferentes combinaciones, elija a o x al azar del conjunto de soluciones válidas para la ecuación anterior. (x debe estar en el rango [l, h] y debe estar (cerca de) un número entero; N * a también debe estar (cerca de) un número entero.

Lior Kogan
fuente