Existe uma maneira eficiente de gerar N inteiros aleatórios em um intervalo que tenha uma determinada soma ou média?

14

Existe uma maneira eficiente de gerar uma combinação aleatória de N números inteiros, como:

  • cada número inteiro está no intervalo [ min, max],
  • os números inteiros têm uma soma de sum,
  • os números inteiros podem aparecer em qualquer ordem (por exemplo, ordem aleatória) e
  • a combinação é escolhida uniformemente aleatoriamente dentre todas as combinações que atendem aos outros requisitos?

Existe um algoritmo semelhante para combinações aleatórias em que os números inteiros devem aparecer na ordem classificada por seus valores (e não em qualquer ordem)?

(Escolher uma combinação apropriada com uma média de meané um caso especial, se sum = N * mean. Esse problema é equivalente a gerar uma partição aleatória uniforme de sumem N partes que estão cada um no intervalo [ min, max] e aparecem em qualquer ordem ou em ordem classificada por seus valores, conforme o caso.)

Estou ciente de que esse problema pode ser resolvido da seguinte maneira para combinações que aparecem em ordem aleatória (EDIT [27 de abril]: Algoritmo modificado.):

  1. Se N * max < sumou N * min > sum, não há solução.

  2. Se N * max == sumhouver apenas uma solução, na qual todos os Nnúmeros são iguais max. Se N * min == sumhouver apenas uma solução, na qual todos os Nnúmeros são iguais min.

  3. Use o algoritmo fornecido em Smith e Tromble ("Sampling from the Unit Simplex", 2004) para gerar N números inteiros aleatórios não negativos com a soma sum - N * min.

  4. Adicione mina cada número gerado dessa maneira.

  5. Se qualquer número for maior que max, vá para a etapa 3.

No entanto, esse algoritmo é lento se maxfor muito menor que sum. Por exemplo, de acordo com meus testes (com uma implementação do caso especial acima envolvendo mean), o algoritmo rejeita, em média,

  • cerca de 1,6 amostras se N = 7, min = 3, max = 10, sum = 42, mas
  • cerca de 30,6 amostras se N = 20, min = 3, max = 10, sum = 120.

Existe uma maneira de modificar esse algoritmo para ser eficiente para N grande, enquanto ainda atende aos requisitos acima?

EDITAR:

Como alternativa sugerida nos comentários, uma maneira eficiente de produzir uma combinação aleatória válida (que satisfaça todos, exceto o último requisito) é:

  1. Calcular X, o número de combinações válidas possíveis sum, dado,, mine max.
  2. Escolha Yum número inteiro aleatório uniforme em [0, X).
  3. Converter ("unrank") Yem uma combinação válida.

No entanto, existe uma fórmula para calcular o número de combinações válidas (ou permutações) e existe uma maneira de converter um número inteiro em uma combinação válida? [EDIT (28 de abril): o mesmo para permutações em vez de combinações].

EDIT (27 de abril):

Depois de ler a Geração de variável aleatória não uniforme de Devroye (1986), posso confirmar que esse é um problema de gerar uma partição aleatória. Além disso, o Exercício 2 (especialmente a parte E) na página 661 é relevante para esta pergunta.

EDIT (28 de abril):

Como se viu, o algoritmo que dei é uniforme, onde os números inteiros envolvidos são dados em ordem aleatória , em oposição à ordem classificada por seus valores . Como ambos os problemas são de interesse geral, modifiquei esta questão para buscar uma resposta canônica para ambos.

O código Ruby a seguir pode ser usado para verificar possíveis soluções de uniformidade (onde algorithm(...)está o 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

EDIT (29 de abril): Re-adicionado código Ruby da implementação atual.

O exemplo de código a seguir é dado em Ruby, mas minha pergunta é independente da linguagem de programação:

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.
fonte
Você poderia esclarecer seu terceiro requisito? Você precisa de uma uniformidade entre todas as combinações possíveis (incluindo aquelas com a média incorreta) ou entre todas as combinações válidas (ou seja, aquelas com a média correta)?
user58697 23/04
Todas as combinações válidas, ou seja, todas as combinações que atendem aos outros requisitos.
Peter O.
Se tivéssemos uma maneira de contar e desagrupar partições de uma soma restrita a N números inteiros em [min, max], a escolha de uma dessas partições aleatoriamente e desagregadas representaria uma distribuição uniforme, e isso seria mais eficiente que o método atual? Qual o tamanho da soma e N?
גלעד ברקן 26/04
Não sei o que você quer dizer com "desmembrar partições de uma soma", e não estou ciente de uma prova de que isso resulta em uma distribuição uniforme no significado desta pergunta. Para esta pergunta, ambos sume Nsão efetivamente ilimitados (dentro da razão). Estou procurando uma resposta canônica porque o problema subjacente aparece em muitas perguntas feitas no Stack Overflow, incluindo esta e esta . @ גלעדברקן
Peter O.
Se atribuirmos a cada combinação possível uma "classificação" (ou índice) em um arranjo ordenado de todas elas, "sem classificação" significaria gerar a combinação, dada sua classificação (e N, min e max, é claro). Por que a escolha de uma dentre todas as combinações possíveis não se conforma a uma distribuição uniforme?
גלעד ברקן 26/04

Respostas:

5

Aqui está a minha solução em Java. É totalmente funcional e contém dois geradores: PermutationPartitionGeneratorpara partições não classificadas e CombinationPartitionGeneratorpara partições classificadas. Seu gerador também foi implementado na classe SmithTromblePartitionGeneratorpara comparação. A classe SequentialEnumeratorenumera todas as partições possíveis (não classificadas ou classificadas, dependendo do parâmetro) em ordem seqüencial. Adicionei testes completos (incluindo seus casos de teste) para todos esses geradores. A implementação é auto-explicável na maior parte. Se você tiver alguma dúvida, eu responderei em alguns dias.

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();
        }
    }
}

Você pode tentar isso no Ideone .

John McClane
fonte
Obrigado pela sua resposta; isso funciona bem. Eu descrevi o gerador de permutação em outra resposta aqui; respondeu outra pergunta com sua ajuda; e em breve incluiremos seu algoritmo no código de exemplo Python do meu artigo sobre métodos de geração aleatória.
Peter O.
Só para ficar claro. Esse algoritmo depende da geração de todas as partições / composições possíveis para fazer a amostra?
Joseph Wood
@ JosephphWood Não, depende de contar todos eles. Isso é feito apenas uma vez na inicialização do gerador e é bastante eficaz porque utiliza a abordagem de programação dinâmica.
John McClane
Como a programação dinâmica pode resolver o problema relacionado de escolher uma partição aleatória uniforme de 'soma' em N números escolhidos aleatoriamente com substituição de uma lista ( exemplo ) ou sem substituição ( exemplo ), ou como esse problema pode ser resolvido?
Peter O.
@PeterO. Você precisa contar todas as partições possíveis pelo mesmo método que no meu algoritmo, mas desta vez você precisa subtrair apenas números permitidos da soma. É muito tempo para comentar, você pode fazer uma pergunta separada. Suspeito que se possa resolver quatro problemas diferentes através da mesma abordagem. Suponha que você tenha uma lista de números inteiros distintos para escolher (este é apenas um intervalo contínuo nesta pergunta). Em seguida, você pode gerar matrizes aleatórias de um determinado comprimento, consistindo em números dessa lista com a soma especificada, se as matrizes forem classificadas / não classificadas e permitir / impedir a repetição.
John McClane
3

Aqui está o algoritmo do PermutationPartitionGenerator de John McClane, em outra resposta nesta página. Possui duas fases, a saber, uma fase de configuração e uma fase de amostragem, e gera nnúmeros aleatórios em [ min, max] com a soma sum, onde os números são listados em ordem aleatória.

Fase de configuração: Primeiro, uma tabela de solução é criada usando as seguintes fórmulas ( t(y, x)onde yestá em [0, n] e xem [0, sum - n * min]):

  • t (0, j) = 1 se j == 0; 0 caso contrário
  • t (i, j) = t (i-1, j) + t (i-1, j-1) + ... + t (i-1, j- (max-min))

Aqui, t (y, x) armazena a probabilidade relativa de que a soma dos ynúmeros (no intervalo apropriado) será igual x. Essa probabilidade é relativa a todos os t (y, x) iguais y.

Fase de amostragem: Aqui geramos uma amostra de nnúmeros. Defina spara sum - n * min, em seguida, para cada posição i, começando comn - 1 e trabalhando para trás em 0:

  • Conjunto v como um número inteiro aleatório em [0, t (i + 1, s)).
  • Defina rcomomin .
  • Subtraia t (i, s) de v .
  • Enquanto vpermanecer 0 ou maior, subtraia t (i, s-1) de v, adicione 1 a re subtraia 1 des .
  • O número na posição ina amostra está definido como r.

EDITAR:

Parece que, com alterações triviais no algoritmo acima, é possível que cada número aleatório use um intervalo separado, em vez de usar o mesmo intervalo para todos eles:

Cada número aleatório nas posições i∈ [0,n ) tem um valor mínimo min (i) e um valor máximo max (i).

Seja adjsum= sum- Σmin (i).

Fase de configuração: Primeiro, uma tabela de solução é criada usando as seguintes fórmulas ( t(y, x)onde yestá em [0, n] e xem [0, adjsum]):

  • t (0, j) = 1 se j == 0; 0 caso contrário
  • t (i, j) = t (i-1, j) + t (i-1, j-1) + ... + t (i-1, j- (max (i-1) -min (i -1)) )

A fase de amostragem é então exatamente a mesma de antes, exceto que configuramos spara adjsum(em vez de sum - n * min) e configuramos rpara min (i) (em vez de min).


EDITAR:

Para CombinationPartitionGenerator, de John McClane, as fases de configuração e amostragem são as seguintes.

Fase de configuração: Primeiro, uma tabela de solução é criada usando as seguintes fórmulas ( t(z, y, x)onde zestá em [0, n], yestá em [0, max - min] e xestá em [0, sum - n * min]):

  • t (0, j, k) = 1 se k == 0; 0 caso contrário
  • 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 amostragem: Aqui geramos uma amostra de nnúmeros. Defina scomo sum - n * mine mrangepara max - min, em seguida, para cada posição i, iniciando n - 1e trabalhando para trás em 0:

  • Defina vcomo um número inteiro aleatório em [0, t (i + 1, intervalo, s)).
  • Definir mrangepara min ( mrange, s)
  • Subtrair mrangede s.
  • Defina rcomomin + mrange .
  • Subtrair t ( i, mrange, s) a partir de v.
  • Enquanto vrestos 0 ou maior, adicionar 1 a s, subtrair 1 re 1 a partir de mrange, em seguida, subtrair t ( i, mrange, s) a partir de v.
  • O número na posição ina amostra está definido como r.
Peter O.
fonte
2

Eu não testei isso, por isso não é realmente uma resposta, apenas algo para tentar que é muito longo para caber em um comentário. Comece com um array que atenda aos dois primeiros critérios e brinque com ele para que ele ainda atenda aos dois primeiros, mas é muito mais aleatório.

Se a média for um número inteiro, sua matriz inicial pode ser [4, 4, 4, ... 4] ou talvez [3, 4, 5, 3, 4, 5, ... 5, 8, 0] ou algo simples assim. Para uma média de 4,5, tente [4, 5, 4, 5, ... 4, 5].

Em seguida, escolha um par de números num1e num2, na matriz. Provavelmente, o primeiro número deve ser tomado em ordem, como no embaralhamento de Fisher-Yates, o segundo número deve ser escolhido aleatoriamente. A ordem do primeiro número garante que todos os números sejam selecionados pelo menos uma vez.

Agora calcule max-num1e num2-min. Essas são as distâncias entre os dois números maxe os minlimites. Defina limitpara a menor das duas distâncias. Essa é a alteração máxima permitida, que não colocará um ou outro número fora dos limites permitidos. Se limitfor zero, pule este par.

Escolha um número inteiro aleatório no intervalo [1, limit]: chame-o change. Eu omito 0 do intervalo selecionável, pois não tem efeito. Os testes podem mostrar que você obtém melhor aleatoriedade ao incluí-la; Não tenho certeza.

Agora defina num1 <- num1 + changee num2 <- num2 - change. Isso não afetará o valor médio e todos os elementos da matriz ainda estão dentro dos limites necessários.

Você precisará percorrer toda a matriz pelo menos uma vez. O teste deve mostrar se você precisa executá-lo mais de uma vez para obter algo suficientemente aleatório.

ETA: incluir 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
fonte
Eu testei e, infelizmente, seu algoritmo não forma uma distribuição uniforme de todas as soluções, não importa quantas iterações eu faça.
Peter O.
Ah bem. Valeu a pena tentar de qualquer maneira. :(
rossum 28/04
1

Como o OP ressalta, a capacidade de desclassificar eficientemente é muito poderosa. Se pudermos fazer isso, a geração de uma distribuição uniforme de partições pode ser feita em três etapas (reafirmando o que o OP estabeleceu na pergunta):

  1. Calcule o número total, M , de partições de comprimento N do número, de summodo que as partes estejam no intervalo [ min, max].
  2. Gere uma distribuição uniforme de números inteiros de [1, M].
  3. Desagrupe cada número inteiro da etapa 2 em sua respectiva partição.

Abaixo, nos concentramos apenas na geração da n- ésima partição, pois há uma quantidade abundante de informações na geração de uma distribuição uniforme de número inteiro em um determinado intervalo. Aqui está um C++algoritmo simples de desagregação que deve ser fácil de traduzir para outros idiomas (NB ainda não descobri como desagregar o caso da composição (por exemplo, a ordem é importante)).

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;
}

A pCountfunção de burro de carga é 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 função é baseada na excelente resposta para Existe um algoritmo eficiente para particionamento inteiro com número restrito de peças? pelo usuário @ m69_snarky_and_unwelcoming. O dado acima é uma ligeira modificação do algoritmo simples (aquele sem memorização). Isso pode ser facilmente modificado para incorporar a memorização para maior eficiência. Por enquanto, deixaremos isso de lado e focaremos na parte sem classificação.

Explicação de unRank

Primeiro, observamos que há um mapeamento individual das partições de comprimento N do número, de summodo que as partes estejam no intervalo [ min, max] até as partições restritas de comprimento N do número sum - N * (min - 1)com partes em [ 1, max - (min - 1)].

Como um pequeno exemplo, considere as partições 50de comprimento 4tais que the min = 10e the max = 15. Isso terá a mesma estrutura que as partições restritas 50 - 4 * (10 - 1) = 14de comprimento 4com a 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

Com isso em mente, para contar facilmente, poderíamos adicionar uma etapa 1a para traduzir o problema para o caso "unit", se desejar.

Agora, simplesmente temos um problema de contagem. Como o @ m69 exibe brilhantemente, a contagem de partições pode ser facilmente obtida dividindo o problema em problemas menores. A função @ m69 fornece nos dá 90% do caminho, só precisamos descobrir o que fazer com a restrição adicional de que existe um limite. É aqui que chegamos:

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

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

Também devemos ter em mente que myMaxisso diminuirá à medida que avançamos. Isso faz sentido se olharmos para a 6 ª partição acima:

2   2   4   6

Para contar o número de partições daqui em diante, devemos continuar aplicando a tradução ao caso "unit". Isso se parece com:

1   1   3   5

Onde, como no passo anterior, tínhamos um máximo de 6, agora consideramos apenas um máximo de 5.

Com isso em mente, desarranjar a partição não é diferente de desarranjar uma permutação ou combinação padrão. Devemos poder contar o número de partições em uma determinada seção. Por exemplo, para contar o número de partições que começam com 10acima, tudo o que fazemos é remover 10a primeira coluna:

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

Traduzir para o caso da unidade:

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

e ligue para pCount:

pCount(13, 3, 6) = 5

Dado um número inteiro aleatório para desagrupar, continuamos calculando o número de partições em seções cada vez menores (como fizemos acima) até preenchermos nosso vetor de índice.

Exemplos

Dada min = 3, max = 10, n = 7, e sum = 42, aqui está um ideone demo que gera 20 partições aleatórias. A saída está abaixo:

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

O índice lexicográfico está à esquerda e a partição não classificada, à direita.

Joseph Wood
fonte
11
Como se vê, essa é uma alternativa muito boa e, de fato, se torna eficiente com a memorização.
Peter O.
0

Se você gerar 0≤a≤1 dos valores aleatórios no intervalo [l, x-1] uniformemente e 1-a dos valores aleatórios no intervalo [x, h] uniformemente, a média esperada seria:

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

Então, se você quer um m específico, pode jogar com a e x.

Por exemplo, se você definir x = m: a = (hm) / (h-l + 1).

Para garantir uma probabilidade mais próxima do uniforme para diferentes combinações, escolha a ou x aleatoriamente do conjunto de soluções válidas para a equação acima. (x deve estar no intervalo [l, h] e deve ser (próximo a) um número inteiro; N * a deve ser (próximo a) um número inteiro também.

Lior Kogan
fonte