vectores usar lenguaje leer imprimir eliminar ejemplos como caracteres cadenas cadena arreglo r performance combinatorics

usar - leer cadena de caracteres en c



Muestra aleatoria de vector de caracteres, sin elementos prefijándose entre sí. (8)

Introducción

Encontré esta pregunta tan interesante que me sentí obligada a reflexionar y, en última instancia, a dar mi propia respuesta. Dado que el algoritmo al que llegué no sigue inmediatamente de la descripción del problema, comenzaré explicando cómo llegué a esta solución y luego proporcionaré una implementación de ejemplo en C ++ (nunca he escrito R).

Desarrollo de la Solución.

A primera vista

La lectura de la descripción del problema fue inicialmente confusa, pero tan pronto como vi la edición con las imágenes de los árboles, comprendí de inmediato el problema y mi intuición sugirió que un árbol binario también era una solución: construir un árbol tamaño 1), y rompa el árbol en una colección de árboles más pequeños después de eliminar las ramas y los antepasados ​​al hacer una selección.

Aunque inicialmente esto parecía bueno, el proceso de selección y el mantenimiento de la colección serían difíciles. Aún así, el árbol parecía que debería jugar un papel importante en cualquier solución.

Revisión 1

No rompas el árbol. En su lugar, tenga una carga útil de datos booleanos en cada nodo que indique si ya se ha eliminado. Esto deja solo un árbol que mantiene la forma.

Sin embargo, tenga en cuenta que esto no es un árbol binario, en realidad es un árbol binario completo de profundidad max_len-1.

Revisión 2

Un árbol binario completo se puede representar muy bien como una matriz. La representación de matriz típica utilizada en la búsqueda en primer lugar del árbol, tiene las siguientes propiedades:

Let x be the array index. x = 0 is the root of the entire tree left_child(x) = 2x + 1 right_child(x) = 2x + 2 parent(x) = floor((n-1)/2)

En la siguiente ilustración, cada nodo está etiquetado con su índice de matriz:

Como una matriz, ocupa menos memoria (no más punteros), la memoria utilizada es contigua (buena para el almacenamiento en caché) y puede ir completamente en la pila en lugar del montón (suponiendo que su idioma le dé una opción). Por supuesto, algunas condiciones se aplican aquí, en particular el tamaño de la matriz. Volveré a eso más tarde.

Al igual que en la Revisión 1, los datos almacenados en la matriz serán booleanos: verdadero para disponible, falso para no disponible. Dado que el nodo raíz no es realmente una opción válida, el índice 0 debe inicializarse en falso. El problema de cómo hacer una selección aún permanece:

A medida que se eliminan los indicios, es trivial hacer un seguimiento de cuántos se eliminan y, por lo tanto, cuántos quedan. Elija un número aleatorio en ese rango, luego recorra la matriz hasta que haya visto que muchos de los indicadores están establecidos en verdadero (se incluye el índice actual). El índice al que llega es entonces la selección a realizar. Seleccione hasta que se seleccionen n indicaciones, o no quede nada para seleccionar.

Este es un algoritmo completo que funcionará, sin embargo, hay espacio para mejorar en el proceso de selección, y también existe un problema práctico de tamaño que no se ha abordado: el tamaño del arreglo sería O (2 ^ n). A medida que n aumenta, primero desaparece el beneficio de almacenamiento en caché, luego los datos comienzan a paginarse en el disco y, en algún momento, resulta imposible almacenarlos.

Revisión 3

Decidí abordar el problema más fácil primero: mejorar el proceso de selección.

Escanear la matriz de izquierda a derecha es un desperdicio. Puede ser más eficiente rastrear los rangos que se han eliminado en lugar de verificar y encontrar varias falsificaciones seguidas. Sin embargo, nuestra representación de árbol no es ideal para esto, ya que pocos de los nodos que serían eliminados en cada ronda serían contiguos en la matriz.

Al reorganizar la forma en que la matriz se asigna al árbol, es posible, sin embargo, explotar mejor esto. En particular, utilicemos la búsqueda de primer orden de profundidad en lugar de la búsqueda de primer alcance. Para hacerlo, el árbol debe ser de tamaño fijo, como es el caso de este problema. También es menos obvio cómo se conectan matemáticamente los índices de los nodos hijo y padre.

Al utilizar este arreglo, se garantiza que cualquier elección que no sea una hoja elimine un rango contiguo: su subárbol.

Revisión 4

Al rastrear los rangos eliminados, ya no hay necesidad de los datos verdaderos / falsos, y por lo tanto, no es necesaria la matriz o el árbol. En cada sorteo aleatorio, los rangos eliminados se pueden usar para encontrar rápidamente el nodo a seleccionar. Todos los antepasados ​​y todo el subárbol se eliminan, y pueden representarse como rangos que pueden combinarse fácilmente con los demás.

La tarea final es convertir los nodos seleccionados en la representación de cadena que el OP deseaba. Eso es bastante fácil, ya que este árbol binario aún mantiene un orden estricto: atravesar desde la raíz, todos los elementos> = el hijo derecho está a la derecha, y los otros están a la izquierda. Por lo tanto, buscar en el árbol proporcionará tanto la lista de ancestros como la cadena binaria al agregar ''0'' cuando se recorra hacia la izquierda; o un ''1'' cuando se cruza a la derecha.

Implementación de la muestra

#include <stdint.h> #include <algorithm> #include <cmath> #include <list> #include <deque> #include <ctime> #include <cstdlib> #include <iostream> /* * A range of values of the form (a, b), where a <= b, and is inclusive. * Ex (1,1) is the range from 1 to 1 (ie: just 1) */ class Range { private: friend bool operator< (const Range& lhs, const Range& rhs); friend std::ostream& operator<<(std::ostream& os, const Range& obj); int64_t m_start; int64_t m_end; public: Range(int64_t start, int64_t end) : m_start(start), m_end(end) {} int64_t getStart() const { return m_start; } int64_t getEnd() const { return m_end; } int64_t size() const { return m_end - m_start + 1; } bool canMerge(const Range& other) const { return !((other.m_start > m_end + 1) || (m_start > other.m_end + 1)); } int64_t merge(const Range& other) { int64_t change = 0; if (m_start > other.m_start) { change += m_start - other.m_start; m_start = other.m_start; } if (other.m_end > m_end) { change += other.m_end - m_end; m_end = other.m_end; } return change; } }; inline bool operator< (const Range& lhs, const Range& rhs){return lhs.m_start < rhs.m_start;} std::ostream& operator<<(std::ostream& os, const Range& obj) { os << ''('' << obj.m_start << '','' << obj.m_end << '')''; return os; } /* * Stuct to allow returning of multiple values */ struct NodeInfo { int64_t subTreeSize; int64_t depth; std::list<int64_t> ancestors; std::string representation; }; /* * Collection of functions representing a complete binary tree * as an array created using pre-order depth-first search, * with 0 as the root. * Depth of the root is defined as 0. */ class Tree { private: int64_t m_depth; public: Tree(int64_t depth) : m_depth(depth) {} int64_t size() const { return (int64_t(1) << (m_depth+1))-1; } int64_t getDepthOf(int64_t node) const{ if (node == 0) { return 0; } int64_t searchDepth = m_depth; int64_t currentDepth = 1; while (true) { int64_t rightChild = int64_t(1) << searchDepth; if (node == 1 || node == rightChild) { break; } else if (node > rightChild) { node -= rightChild; } else { node -= 1; } currentDepth += 1; searchDepth -= 1; } return currentDepth; } int64_t getSubtreeSizeOf(int64_t node, int64_t nodeDepth = -1) const { if (node == 0) { return size(); } if (nodeDepth == -1) { nodeDepth = getDepthOf(node); } return (int64_t(1) << (m_depth + 1 - nodeDepth)) - 1; } int64_t getLeftChildOf(int64_t node, int64_t nodeDepth = -1) const { if (nodeDepth == -1) { nodeDepth = getDepthOf(node); } if (nodeDepth == m_depth) { return -1; } return node + 1; } int64_t getRightChildOf(int64_t node, int64_t nodeDepth = -1) const { if (nodeDepth == -1) { nodeDepth = getDepthOf(node); } if (nodeDepth == m_depth) { return -1; } return node + 1 + ((getSubtreeSizeOf(node, nodeDepth) - 1) / 2); } NodeInfo getNodeInfo(int64_t node) const { NodeInfo info; int64_t depth = 0; int64_t currentNode = 0; while (currentNode != node) { if (currentNode != 0) { info.ancestors.push_back(currentNode); } int64_t rightChild = getRightChildOf(currentNode, depth); if (rightChild == -1) { break; } else if (node >= rightChild) { info.representation += ''1''; currentNode = rightChild; } else { info.representation += ''0''; currentNode = getLeftChildOf(currentNode, depth); } depth++; } info.depth = depth; info.subTreeSize = getSubtreeSizeOf(node, depth); return info; } }; // random selection amongst remaining allowed nodes int64_t selectNode(const std::deque<Range>& eliminationList, int64_t poolSize, std::mt19937_64& randomGenerator) { std::uniform_int_distribution<> randomDistribution(1, poolSize); int64_t selection = randomDistribution(randomGenerator); for (auto const& range : eliminationList) { if (selection >= range.getStart()) { selection += range.size(); } else { break; } } return selection; } // determin how many nodes have been elimintated int64_t countEliminated(const std::deque<Range>& eliminationList) { int64_t count = 0; for (auto const& range : eliminationList) { count += range.size(); } return count; } // merge all the elimination ranges to listA, and return the number of new elimintations int64_t mergeEliminations(std::deque<Range>& listA, std::deque<Range>& listB) { if(listB.empty()) { return 0; } if(listA.empty()) { listA.swap(listB); return countEliminated(listA); } int64_t newEliminations = 0; int64_t x = 0; auto listA_iter = listA.begin(); auto listB_iter = listB.begin(); while (listB_iter != listB.end()) { if (listA_iter == listA.end()) { listA_iter = listA.insert(listA_iter, *listB_iter); x = listB_iter->size(); assert(x >= 0); newEliminations += x; ++listB_iter; } else if (listA_iter->canMerge(*listB_iter)) { x = listA_iter->merge(*listB_iter); assert(x >= 0); newEliminations += x; ++listB_iter; } else if (*listB_iter < *listA_iter) { listA_iter = listA.insert(listA_iter, *listB_iter) + 1; x = listB_iter->size(); assert(x >= 0); newEliminations += x; ++listB_iter; } else if ((listA_iter+1) != listA.end() && listA_iter->canMerge(*(listA_iter+1))) { listA_iter->merge(*(listA_iter+1)); listA_iter = listA.erase(listA_iter+1); } else { ++listA_iter; } } while (listA_iter != listA.end()) { if ((listA_iter+1) != listA.end() && listA_iter->canMerge(*(listA_iter+1))) { listA_iter->merge(*(listA_iter+1)); listA_iter = listA.erase(listA_iter+1); } else { ++listA_iter; } } return newEliminations; } int main (int argc, char** argv) { std::random_device rd; std::mt19937_64 randomGenerator(rd()); int64_t max_len = std::stoll(argv[1]); int64_t num_samples = std::stoll(argv[2]); int64_t samplesRemaining = num_samples; Tree tree(max_len); int64_t poolSize = tree.size() - 1; std::deque<Range> eliminationList; std::deque<Range> eliminated; std::list<std::string> foundList; while (samplesRemaining > 0 && poolSize > 0) { // find a valid node int64_t selectedNode = selectNode(eliminationList, poolSize, randomGenerator); NodeInfo info = tree.getNodeInfo(selectedNode); foundList.push_back(info.representation); samplesRemaining--; // determine which nodes this choice eliminates eliminated.clear(); for( auto const& ancestor : info.ancestors) { Range r(ancestor, ancestor); if(eliminated.empty() || !eliminated.back().canMerge(r)) { eliminated.push_back(r); } else { eliminated.back().merge(r); } } Range r(selectedNode, selectedNode + info.subTreeSize - 1); if(eliminated.empty() || !eliminated.back().canMerge(r)) { eliminated.push_back(r); } else { eliminated.back().merge(r); } // add the eliminated nodes to the existing list poolSize -= mergeEliminations(eliminationList, eliminated); } // Print some stats // std::cout << "tree: " << tree.size() << " samplesRemaining: " // << samplesRemaining << " poolSize: " // << poolSize << " samples: " << foundList.size() // << " eliminated: " // << countEliminated(eliminationList) << std::endl; // Print list of binary strings // std::cout << "list:"; // for (auto const& s : foundList) { // std::cout << " " << s; // } // std::cout << std::endl; }

Pensamientos adicionales

Este algoritmo se escalará extremadamente bien para max_len. La escala con n no es muy buena, pero según mi propio perfil parece funcionar mejor que las otras soluciones.

Este algoritmo podría modificarse con poco esfuerzo para permitir cadenas que contengan más que solo ''0'' y ''1''. Más letras posibles en las palabras aumentan el abanico de salida del árbol y eliminarían un rango más amplio con cada selección, aún con todos los nodos en cada subárbol que permanecen contiguos.

Considere un vector de caracteres, pool , cuyos elementos son números binarios (con relleno cero) con hasta max_len dígitos.

max_len <- 4 pool <- unlist(lapply(seq_len(max_len), function(x) do.call(paste0, expand.grid(rep(list(c(''0'', ''1'')), x))))) pool ## [1] "0" "1" "00" "10" "01" "11" "000" "100" "010" "110" ## [11] "001" "101" "011" "111" "0000" "1000" "0100" "1100" "0010" "1010" ## [21] "0110" "1110" "0001" "1001" "0101" "1101" "0011" "1011" "0111" "1111"

Me gustaría muestrear n de estos elementos, con la restricción de que ninguno de los elementos muestreados son prefixes de ninguno de los otros elementos muestreados (es decir, si muestreamos 1101 , prohibimos 1 , 11 y 110 , mientras que si muestreamos 1 , prohibimos aquellos elementos que comienzan con 1 , como 10 , 11 , 100 , etc.).

El siguiente es mi intento de usar while , pero por supuesto esto es lento cuando n es grande (o cuando se acerca a 2^max_len ).

set.seed(1) n <- 10 chosen <- sample(pool, n) while(any(rowSums(outer(paste0(''^'', chosen), chosen, Vectorize(grepl))) > 1)) { prefixes <- rowSums(outer(paste0(''^'', chosen), chosen, Vectorize(grepl))) > 1 pool <- pool[rowSums(Vectorize(grepl, ''pattern'')( paste0(''^'', chosen[!prefixes]), pool)) == 0] chosen <- c(chosen[!prefixes], sample(pool, sum(prefixes))) } chosen ## [1] "0100" "0101" "0001" "0011" "1000" "111" "0000" "0110" "1100" "0111"

Esto se puede mejorar ligeramente eliminando inicialmente del pool aquellos elementos cuya inclusión significaría que no hay elementos suficientes en el pool para tomar una muestra total de tamaño n . Por ejemplo, cuando max_len = 4 y n > 9 , podemos eliminar inmediatamente 0 y 1 de la pool , ya que al incluir cualquiera de ellos, la muestra máxima sería 9 (ya sea 0 y los ocho elementos de 4 caracteres que comienzan con 1 , o 1 y Ocho elementos de 4 caracteres que comienzan con 0 ).

Según esta lógica, podríamos omitir elementos de la pool , antes de tomar la muestra inicial, con algo como:

pool <- pool[ nchar(pool) > tail(which(n > (2^max_len - rev(2^(0:max_len))[-1] + 1)), 1)]

¿Alguien puede pensar en un mejor enfoque? Siento que estoy pasando por alto algo mucho más simple.

EDITAR

Para aclarar mis intenciones, retrataré el conjunto como un conjunto de ramas, donde las uniones y las puntas son los nodos (elementos del pool ). Supongamos que se dibujó el nodo amarillo en la siguiente figura (es decir, 010). Ahora, toda la "rama" roja, que consta de los nodos 0, 01 y 010, se elimina de la agrupación. Esto es lo que quise decir al prohibir el muestreo de nodos que "prefijan" los nodos que ya están en nuestra muestra (así como los que ya están prefijados por los nodos en nuestra muestra).

Si el nodo muestreado estaba a mitad de camino a lo largo de una rama, como 01 en la siguiente figura, entonces todos los nodos rojos (0, 01, 010 y 011) no están permitidos, ya que 0 prefijos 01 y 01 prefijos tanto 010 como 011.

No quiero probar 1 o 0 en cada cruce (es decir, caminar a lo largo de las ramas lanzando monedas en los tenedores); está bien tener ambos en la muestra, siempre y cuando: (1) los padres (o los abuelos, etc.) .) o los hijos (nietos, etc.) del nodo no se han muestreado; y (2) al muestrear el nodo, habrá suficientes nodos restantes para lograr la muestra deseada de tamaño n .

En la segunda figura anterior, si 010 fue la primera selección, todos los nodos en los nodos negros siguen siendo (actualmente) válidos, suponiendo que n <= 4 . Por ejemplo, si n==4 y muestreamos el nodo 1 a continuación (y por lo tanto, nuestras selecciones ahora incluían 01 y 1), posteriormente rechazaremos el nodo 00 (debido a la regla 2 anterior), pero aún podríamos seleccionar 000 y 001, lo que nos da nuestra Muestra de 4 elementos. Si n==5 , por otro lado, el nodo 1 se rechazaría en esta etapa.


Introducción

Esta es una variación numérica del algoritmo de cadena que implementamos en la otra respuesta. Es más rápido y no requiere crear ni ordenar la agrupación.

Esquema del algoritmo

Podemos usar números enteros para representar sus cadenas binarias, lo que simplifica enormemente el problema de la generación de agrupaciones y la eliminación secuencial de valores. Por ejemplo, con max_len==3 , podemos tomar el número 1-- (donde - representa el relleno) para significar 4 en decimal. Además, podemos determinar que los números que requieren eliminación si seleccionamos este número son aquellos entre 4 y 4 + 2 ^ x - 1 . Aquí x es el número de elementos de relleno (2 en este caso), por lo que los números a eliminar están entre 4 y 4 + 2 ^ 2 - 1 (o entre 4 y 7 , expresados ​​como 100 , 110 y 111 ).

Para poder resolver exactamente su problema, necesitamos un poco de arrugas, ya que los números que posiblemente serían los mismos en binario podrían ser distintos para algunas partes de su algoritmo. Por ejemplo, 100 , 10- y 1-- son todos el mismo número, pero deben tratarse de manera diferente en su esquema. En un mundo max_len==3 , tenemos 8 números posibles, pero 14 representaciones posibles:

0 - 000: 0--, 00- 1 - 001: 2 - 010: 01- 3 - 011: 4 - 100: 1--, 10- 5 - 101: 6 - 110: 11- 7 - 111:

Así que 0 y 4 tienen tres codificaciones posibles, 2 y 6 tienen dos, y todas las demás son solo una. Necesitamos generar un grupo de enteros que represente la mayor probabilidad de selección para los números con representación múltiple, así como un mecanismo para rastrear cuántos espacios en blanco incluye el número. Podemos hacer esto agregando algunos bits al final del número para indicar las ponderaciones que queremos. Así que nuestros números se vuelven (usamos dos bits aquí):

jbaum | int | bin | bin.enc | int.enc 0-- | 0 | 000 | 00000 | 0 00- | 0 | 000 | 00001 | 1 000 | 0 | 000 | 00010 | 2 001 | 1 | 001 | 00100 | 3 01- | 2 | 010 | 01000 | 4 010 | 2 | 010 | 01001 | 5 011 | 3 | 011 | 01101 | 6 1-- | 4 | 100 | 10000 | 7 10- | 4 | 100 | 10001 | 8 100 | 4 | 100 | 10010 | 9 101 | 5 | 101 | 10100 | 10 11- | 6 | 110 | 11000 | 11 110 | 6 | 110 | 11001 | 12 111 | 7 | 111 | 11100 | 13

Algunas propiedades útiles:

  • enc.bits representa la cantidad de bits que necesitamos para la codificación (dos en este caso)
  • int.enc %% enc.bits nos dice cuántos de los números están explícitamente especificados
  • int.enc %/% enc.bits devuelve int
  • int * 2 ^ enc.bits + explicitly.specified devuelve int.enc

Tenga en cuenta que explicitly.specified especifica aquí entre 0 y max_len - 1 en nuestra implementación, ya que siempre hay al menos un dígito especificado. Ahora tenemos una codificación que representa completamente su estructura de datos utilizando solo números enteros. Podemos muestrear desde enteros y reproducir los resultados deseados, con los pesos correctos, etc. Una limitación de este enfoque es que usamos enteros de 32 bits en R, y debemos reservar algunos bits para la codificación, por lo que nos limitamos a grupos de max_len==25 o menos. Usted podría ir más grande si usó números enteros especificados por puntos flotantes de doble precisión, pero no lo hicimos aquí.

Evitando Selecciones Duplicadas

Hay dos formas básicas para asegurarnos de que no elegimos el mismo valor dos veces

  1. Haga un seguimiento de los valores que quedan disponibles para seleccionar y muestre de forma aleatoria entre ellos
  2. Muestra aleatoriamente de todos los valores posibles, y luego verifica si el valor ha sido seleccionado previamente, y si lo ha hecho, muestra nuevamente

Si bien la primera opción parece la más limpia, en realidad es muy costosa computacionalmente. Se requiere realizar una exploración vectorial de todos los valores posibles para cada selección para descalificar previamente los valores seleccionados o crear un vector de reducción que contenga valores no descalificados. La opción de reducción solo es más eficiente que el escaneo vectorial si se hace que el vector se reduzca por referencia a través del código C, pero incluso entonces requiere traducciones repetidas de porciones potencialmente grandes del vector, y requiere C.

Aquí usamos el método # 2. Esto nos permite mezclar aleatoriamente el universo de valores posibles una vez, y luego seleccionar cada valor en secuencia, verificar que no haya sido descalificado, y si lo ha hecho, elegir otro, etc. Esto funciona porque es trivial verificar si un el valor ha sido elegido como resultado de nuestra codificación de valor; podemos inferir la ubicación de un valor en una tabla ordenada en base únicamente al valor . Por lo tanto, registramos el estado de cada valor en una tabla ordenada, y podemos actualizar o buscar ese estado a través del acceso directo al índice (no se requiere escaneo).

Ejemplos

La implementación de este algoritmo en base R está disponible como una esencia . Esta implementación particular solo extrae los sorteos completos. Aquí hay una muestra de 10 dibujos de 8 elementos de un grupo max_len==4 :

# each column represents a draw from a `max_len==4` pool set.seed(6); replicate(10, sample0110b(4, 8)) [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [1,] "1000" "1" "0011" "0010" "100" "0011" "0" "011" "0100" "1011" [2,] "111" "0000" "1101" "0000" "0110" "0100" "1000" "00" "0101" "1001" [3,] "0011" "0110" "1001" "0100" "0000" "0101" "1101" "1111" "10" "1100" [4,] "0100" "0010" "0000" "0101" "1101" "101" "1011" "1101" "0110" "1101" [5,] "101" "0100" "1100" "1100" "0101" "1001" "1001" "1000" "1111" "1111" [6,] "110" "0111" "1011" "111" "1011" "110" "1111" "0100" "0011" "000" [7,] "0101" "0101" "111" "011" "1010" "1000" "1100" "101" "0001" "0101" [8,] "011" "0001" "01" "1010" "0011" "1110" "1110" "1001" "110" "1000"

Originalmente teníamos dos implementaciones que dependían del método # 1 para evitar duplicados, una en la base R y otra en C, pero incluso la versión C no es tan rápida como la nueva versión de la base R cuando n es grande. Estas funciones implementan la capacidad de dibujar sorteos incompletos, por lo que les ofrecemos aquí como referencia:

Puntos de referencia comparativos

Aquí hay un conjunto de puntos de referencia que comparan varias de las funciones que se muestran en esta Q / A. Tiempos en milisegundos. La versión brodie.b es la que se describe en esta respuesta. brodie es la implementación original, brodie.C es el original con algo de C. Todos estos hacen cumplir el requisito de muestras completas. brodie.str es la versión basada en cadena en la otra respuesta.

size n jbaum josilber frank tensibai brodie.b brodie brodie.C brodie.str 1 4 10 11 1 3 1 1 1 1 0 2 4 50 - - - 1 - - - 1 3 4 100 - - - 1 - - - 0 4 4 256 - - - 1 - - - 1 5 4 1000 - - - 1 - - - 1 6 8 10 1 290 6 3 2 2 1 1 7 8 50 388 - 8 8 3 4 3 4 8 8 100 2,506 - 13 18 6 7 5 5 9 8 256 - - 22 27 13 14 12 6 10 8 1000 - - - 27 - - - 7 11 16 10 - - 615 688 31 61 19 424 12 16 50 - - 2,123 2,497 28 276 19 1,764 13 16 100 - - 4,202 4,807 30 451 23 3,166 14 16 256 - - 11,822 11,942 40 1,077 43 8,717 15 16 1000 - - 38,132 44,591 83 3,345 130 27,768

Esto se adapta relativamente bien a piscinas más grandes.

system.time(sample0110b(18, 100000)) user system elapsed 8.441 0.079 8.527

Notas de referencia:

  • frank, y brodie (menos brodie.str) no requieren ninguna agrupación previa, lo que afectará las comparaciones (ver más abajo)
  • Josilber es la versión LP
  • jbaum es el ejemplo OP
  • tensibai se modifica ligeramente para salir en lugar de fallar si el grupo está vacío
  • no está configurado para ejecutar Python, por lo que no puedo comparar / tener en cuenta el almacenamiento en búfer
  • - Representa opciones inviables o demasiado lentas para calcular el tiempo razonablemente.

Los tiempos no incluyen dibujar las agrupaciones ( 0.8 , 2.5 , 401 milisegundos respectivamente para el tamaño 4 , 8 y 16 ), lo que es necesario para las jbaum , josilber y brodie.str , ni clasificarlas ( 0.1 , 2.7 , 3700 milisegundos para el tamaño 4 , 8 y 16 ), que es necesario para brodie.str además del sorteo. Si desea incluirlos o no, depende de cuántas veces ejecute la función para un grupo en particular. Además, es casi seguro que hay mejores formas de generar / clasificar las agrupaciones.

Estos son el tiempo medio de tres carreras con microbenchmark . El código está disponible como una esencia , aunque tenga en cuenta que deberá cargar las sample0110b , sample0110 , sample01101 y sample01 antemano.


Encontré el problema interesante, así que intenté y obtuve esto con habilidades realmente bajas en R (para que esto pueda mejorarse):

Frank versión editada, gracias al consejo de Frank :

library(microbenchmark) library(lineprof) max_len <- 16 pool <- unlist(lapply(seq_len(max_len), function(x) do.call(paste0, expand.grid(rep(list(c(''0'', ''1'')), x))))) n<-100 library(stringr) tree_sample <- function(samples,pool) { results <- vector("integer",samples) # Will be used on a regular basis, compute it in advance PoolLen <- str_length(pool) # Make a mask vector based on the length of each entry of the pool masks <- strtoi(str_pad(str_pad("1",PoolLen,"right","1"),max_len,"right","0"),base=2) # Make an integer vector from "0" right padded orignal: for max_len=4 and pool entry "1" we get "1000" => 8 # This will allow to find this entry as parent of 10 and 11 which become "1000" and "1100", as integer 8 and 12 respectively # once bitwise "anded" with the repective mask "1000" the first bit is striclty the same, so it''s a parent. integerPool <- strtoi(str_pad(pool,max_len,"right","0"),base=2) # Create a vector to filter the available value to sample ok <- rep(TRUE,length(pool)) #Precompute the result of the bitwise and betwwen our integer pool and the masks MaskedPool <- bitwAnd(integerPool,masks) while(samples) { samp <- sample(pool[ok],1) # Get a sample results[samples] <- samp # Store it as result ok[pool == samp] <- FALSE # Remove it from available entries vsamp <- strtoi(str_pad(samp,max_len,"right","0"),base=2) # Get the integer value of the "0" right padded sample mlen <- str_length(samp) # Get sample len #Creation of unitary mask to remove childs of sample mask <- strtoi(paste0(rep(1:0,c(mlen,max_len-mlen)),collapse=""),base=2) # Get the result of bitwise And between the integerPool and the sample mask FilterVec <- bitwAnd(integerPool,mask) # Get the bitwise and result of the sample and it''s mask Childm <- bitwAnd(vsamp,mask) ok[FilterVec == Childm] <- FALSE # Remove from available entries the childs of the sample ok[MaskedPool == bitwAnd(vsamp,masks)] <- FALSE # compare the sample with all the masks to remove parents matching samples <- samples -1 } print(results) } microbenchmark(tree_sample(n,pool),times=10L)

La idea principal es utilizar la comparación de máscara de bits para saber si una muestra es padre de la otra (parte de bit común), si es así, suprima este elemento de la agrupación.

Ahora se necesitan 1,4 s para extraer 100 muestras de un grupo de longitud 16 en mi máquina.


Puede ordenar el grupo para ayudar a decidir qué elementos descalificar. Por ejemplo, mirando un grupo ordenado de tres elementos:

[1] "0" "00" "000" "001" "01" "010" "011" "1" "10" "100" "101" "11" [13] "110" "111"

Puedo decir que puedo descalificar cualquier cosa que siga a mi elemento seleccionado que tenga más caracteres que mi elemento hasta el primer elemento que tenga el mismo número o menos caracteres. Por ejemplo, si selecciono "01", puedo ver inmediatamente que los siguientes dos elementos ("010", "011") deben eliminarse, pero no el siguiente después de eso porque "1" tiene menos caracteres. Eliminar el "0" después es fácil. Aquí hay una implementación:

library(fastmatch) # could use `match`, but we repeatedly search against same hash # `pool` must be sorted! sample01 <- function(pool, n) { picked <- logical(length(pool)) chrs <- nchar(pool) pick.list <- character(n) pool.seq <- seq_along(pool) for(i in seq(n)) { # Make sure pool not exhausted left <- which(!picked) left.len <- length(left) if(!length(left)) break # Sample from pool seq.left <- seq.int(left) pool.left <- pool[left] chrs.left <- chrs[left] pick <- sample(length(pool.left), 1L) # Find all the elements with more characters that are disqualified # and store their indices in `valid` (bad name...) valid.tmp <- chrs.left > chrs.left[[pick]] & seq.left > pick first.invalid <- which(!valid.tmp & seq.left > pick) valid <- if(length(first.invalid)) { pick:(first.invalid[[1L]] - 1L) } else pick:left.len # Translate back to original pool indices since we''re working on a # subset in `pool.left` pool.seq.left <- pool.seq[left] pool.idx <- pool.seq.left[valid] val <- pool[[pool.idx[[1L]]]] # Record the picked value, and all the disqualifications pick.list[[i]] <- val picked[pool.idx] <- TRUE # Disqualify shorter matches to.rem <- vapply( seq.int(nchar(val) - 1), substr, character(1L), x=val, start=1L ) to.rem.idx <- fmatch(to.rem, pool, nomatch=0) picked[to.rem.idx] <- TRUE } pick.list }

Y una función para hacer grupos ordenados (exactamente como su código, pero los resultados ordenados):

make_pool <- function(size) sort( unlist( lapply( seq_len(size), function(x) do.call(paste0, expand.grid(rep(list(c(''0'', ''1'')), x))) ) ) )

Luego, utilizando el grupo max_len 3 (útil para inspeccionar visualmente las cosas se comportan como se espera):

pool3 <- make_pool(3) set.seed(1) sample01(pool3, 8) # [1] "001" "1" "010" "011" "000" "" "" "" sample01(pool3, 8) # [1] "110" "111" "011" "10" "00" "" "" "" sample01(pool3, 8) # [1] "000" "01" "11" "10" "001" "" "" "" sample01(pool3, 8) # [1] "011" "101" "111" "001" "110" "100" "000" "010"

Observe que en el último caso obtenemos todas las combinaciones binarias de 3 dígitos (2 ^ 3) porque por casualidad seguimos muestreando de los de 3 dígitos. Además, con solo un grupo de 3 tamaños hay muchos muestreos que impiden un total de 8 sorteos; usted podría abordar esto con su sugerencia de eliminar las combinaciones que impiden el sorteo completo de la piscina.

Esto es bastante rápido. Mirando el ejemplo de max_len==9 que tomó 2 segundos con la solución alternativa:

pool9 <- make_pool(9) microbenchmark(sample01(pool9, 4)) # Unit: microseconds # expr min lq median uq max neval # sample01(pool9, 4) 493.107 565.015 571.624 593.791 983.663 100

Alrededor de medio milisegundo. También puedes intentar razonablemente piscinas bastante grandes:

pool16 <- make_pool(16) # 131K entries system.time(sample01(pool16, 100)) # user system elapsed # 3.407 0.146 3.552

Esto no es muy rápido, pero estamos hablando de un grupo con elementos de 130K. También hay potencial para la optimización adicional.

Tenga en cuenta que el paso de clasificación es relativamente lento para los grupos grandes, pero no lo estoy contando porque solo necesita hacerlo una vez, y probablemente pueda encontrar un algoritmo razonable para producir el grupo ordenado previamente.

También existe la posibilidad de un enfoque mucho más rápido de enteros a binarios que exploré en una respuesta ahora eliminada, pero eso requiere un poco más de trabajo para vincular exactamente lo que está buscando.


Está en python y no en r, pero jbaums dijo que estaría bien.

Así que aquí está mi contribución, vea los comentarios en la fuente para obtener una explicación de las partes cruciales.
Todavía estoy trabajando en la solución analítica para determinar la cantidad de posibles combinaciones c para un árbol de muestras t y S de profundidad, por lo que puedo mejorar los combs la función. Tal vez alguien lo tiene? Este es realmente el cuello de botella ahora.

El muestreo de 100 nodos de un árbol con profundidad 16 toma alrededor de 8 ms en mi computadora portátil. No es la primera vez, pero se acelera hasta cierto punto cuanto más se muestrea porque se llena combBuffer.

import random class Tree(object): """ :param level: The distance of this node from the root. :type level: int :param parent: This trees parent node :type parent: Tree :param isleft: Determines if this is a left or a right child node. Can be omitted if this is the root node. :type isleft: bool A binary tree representing possible strings which match r''[01]{1,n}''. Its purpose is to be able to sample n of its nodes where none of the sampled nodes'' ids is a prefix for another one. It is possible to change Tree.maxdepth and then reuse the root. All children are created ON DEMAND, which means everything is lazily evaluated. If the Tree gets too big anyway, you can call ''prune'' on any node to delete its children. >>> t = Tree() >>> t.sample(8, toString=True, depth=3) [''111'', ''110'', ''101'', ''100'', ''011'', ''010'', ''001'', ''000''] >>> Tree.maxdepth = 2 >>> t.sample(4, toString=True) [''11'', ''10'', ''01'', ''00''] """ maxdepth = 10 _combBuffer = {} def __init__(self, level=0, parent=None, isleft=None): self.parent = parent self.level = level self.isleft = isleft self._left = None self._right = None @classmethod def setMaxdepth(cls, depth): """ :param depth: The new depth :type depth: int Sets the maxdepth of the Trees. This basically is the depth of the root node. """ if cls.maxdepth == depth: return cls.maxdepth = depth @property def left(self): """This tree''s left child, ''None'' if this is a leave node""" if self.depth == 0: return None if self._left is None: self._left = Tree(self.level+1, self, True) return self._left @property def right(self): """This tree''s right child, ''None'' if this is a leave node""" if self.depth == 0: return None if self._right is None: self._right = Tree(self.level+1, self, False) return self._right @property def depth(self): """ This tree''s depth. (maxdepth-level) """ return self.maxdepth-self.level @property def id(self): """ This tree''s id, string of ''0''s and ''1''s equal to the path from the root to this subtree. Where ''1'' means going left and ''0'' means going right. """ # level 0 is the root node, it has no id if self.level == 0: return '''' # This takes at most Tree.maxdepth recursions. Therefore # it is save to do it this way. We could also save each nodes # id once it is created to avoid recreating it every time, however # this won''t save much time but use quite some space. return self.parent.id + (''1'' if self.isleft else ''0'') @property def leaves(self): """ The amount of leave nodes, this tree has. (2**depth) """ return 2**self.depth def __str__(self): return self.id def __len__(self): return 2*self.leaves-1 def prune(self): """ Recursively prune this tree''s children. """ if self._left is not None: self._left.prune() self._left.parent = None self._left = None if self._right is not None: self._right.prune() self._right.parent = None self._right = None def combs(self, n): """ :param n: The amount of samples to be taken from this tree :type n: int :returns: The amount of possible combinations to choose n samples from this tree Determines recursively the amount of combinations of n nodes to be sampled from this tree. Subsequent calls with same n on trees with same depth will return the result from the previous computation rather than computing it again. >>> t = Tree() >>> Tree.maxdepth = 4 >>> t.combs(16) 1 >>> Tree.maxdepth = 3 >>> t.combs(6) 58 """ # important for the amount of combinations is only n and the depth of # this tree key = (self.depth, n) # We use the dict to save computation time. Calling the function with # equal values on equal nodes just returns the alrady computed value if # possible. if key not in Tree._combBuffer: leaves = self.leaves if n < 0: N = 0 elif n == 0 or self.depth == 0 or n == leaves: N = 1 elif n == 1: return (2*leaves-1) else: if n > leaves/2: # if n > leaves/2, at least n-leaves/2 have to stay on # either side, otherweise the other one would have to # sample more nodes than possible. nMin = n-leaves/2 else: nMin = 0 # The rest n-2*nMin is the amount of samples that are free to # fall on either side free = n-2*nMin N = 0 # sum up the combinations of all possible splits for addLeft in range(0, free+1): nLeft = nMin + addLeft nRight = n - nLeft N += self.left.combs(nLeft)*self.right.combs(nRight) Tree._combBuffer[key] = N return N return Tree._combBuffer[key] def sample(self, n, toString=False, depth=None): """ :param n: How may samples to take from this tree :type n: int :param toString: If ''True'' result will direclty be turned into a list of strings :type toString: bool :param depth: If not None, will overwrite Tree.maxdepth :type depth: int or None :returns: List of n nodes sampled from this tree :throws ValueError: when n is invalid Takes n random samples from this tree where none of the sample''s ids is a prefix for another one''s. For an example see Tree''s docstring. """ if depth is not None: Tree.setMaxdepth(depth) if toString: return [str(e) for e in self.sample(n)] if n < 0: raise ValueError(''Negative sample size is not possible!'') if n == 0: return [] leaves = self.leaves if n > leaves: raise ValueError((''Cannot sample {} nodes, with only {} '' + ''leaves!'').format(n, leaves)) # Only one sample to choose, that is nice! We are free to take any node # from this tree, including this very node itself. if n == 1 and self.level > 0: # This tree has 2*leaves-1 nodes, therefore # the probability that we keep the root node has to be # 1/(2*leaves-1) = P_root. Lets create a random number from the # interval [0, 2*leaves-1). # It will be 0 with probability 1/(2*leaves-1) P_root = random.randint(0, len(self)-1) if P_root == 0: return [self] else: # The probability to land here is 1-P_root # A child tree''s size is (leaves-1) and since it obeys the same # rule as above, the probability for each of its nodes to # ''survive'' is 1/(leaves-1) = P_child. # However all nodes must have equal probability, therefore to # make sure that their probability is also P_root we multiply # them by 1/2*(1-P_root). The latter is already done, the # former will be achieved by the next condition. # If we do everything right, this should hold: # 1/2 * (1-P_root) * P_child = P_root # Lets see... # 1/2 * (1-1/(2*leaves-1)) * (1/leaves-1) # (1-1/(2*leaves-1)) * (1/(2*(leaves-1))) # (1-1/(2*leaves-1)) * (1/(2*leaves-2)) # (1/(2*leaves-2)) - 1/((2*leaves-2) * (2*leaves-1)) # (2*leaves-1)/((2*leaves-2) * (2*leaves-1)) - 1/((2*leaves-2) * (2*leaves-1)) # (2*leaves-2)/((2*leaves-2) * (2*leaves-1)) # 1/(2*leaves-1) # There we go! if random.random() < 0.5: return self.right.sample(1) else: return self.left.sample(1) # Now comes the tricky part... n > 1 therefore we are NOT going to # sample this node. Its probability to be chosen is 0! # It HAS to be 0 since we are definitely sampling from one of its # children which means that this node will be blocked by those samples. # The difficult part now is to prove that the sampling the way we do it # is really random. if n > leaves/2: # if n > leaves/2, at least n-leaves/2 have to stay on either # side, otherweise the other one would have to sample more # nodes than possible. nMin = n-leaves/2 else: nMin = 0 # The rest n-2*nMin is the amount of samples that are free to fall # on either side free = n-2*nMin # Let''s have a look at an example, suppose we were to distribute 5 # samples among two children which have 4 leaves each. # Each child has to get at least 1 sample, so the free samples are 3. # There are 4 different ways to split the samples among the # children (left, right): # (1, 4), (2, 3), (3, 2), (4, 1) # The amount of unique sample combinations per child are # (7, 1), (11, 6), (6, 11), (1, 7) # The amount of total unique samples per possible split are # 7 , 66 , 66 , 7 # In case of the first and last split, all samples have a probability # of 1/7, this was already proven above. # Lets suppose we are good to go and the per sample probabilities for # the other two cases are (1/11, 1/6) and (1/6, 1/11), this way the # overall per sample probabilities for the splits would be: # 1/7 , 1/66 , 1/66 , 1/7 # If we used uniform random to determine the split, all splits would be # equally probable and therefore be multiplied with the same value (1/4) # But this would mean that NOT every sample is equally probable! # We need to know in advance how many sample combinations there will be # for a given split in order to find out the probability to choose it. # In fact, due to the restrictions, this becomes very nasty to # determine. So instead of solving it analytically, I do it numerically # with the method ''combs''. It gives me the amount of possible sample # combinations for a certain amount of samples and a given tree depth. # It will return 146 for this node and 7 for the outer and 66 for the # inner splits. # What we now do is, we take a number from [0, 146). # if it is smaller than 7, we sample from the first split, # if it is smaller than 7+66, we sample from the second split, # ... # This way we get the probabilities we need. r = random.randint(0, self.combs(n)-1) p = 0 for addLeft in xrange(0, free+1): nLeft = nMin + addLeft nRight = n - nLeft p += (self.left.combs(nLeft) * self.right.combs(nRight)) if r < p: return self.left.sample(nLeft) + self.right.sample(nRight) assert False, (''Something really strange happend, p did not sum up '' + ''to combs or r was too big'') def main(): """ Do a microbenchmark. """ import timeit i = 1 main.t = Tree() template = '' {:>2} {:>5} {:>4} {:<5}'' print(template.format(''i'', ''depth'', ''n'', ''time (ms)'')) N = 100 for depth in [4, 8, 15, 16, 17, 18]: for n in [10, 50, 100, 150]: if n > 2**depth: time = ''--'' else: time = timeit.timeit( ''main.t.sample({}, depth={})''.format(n, depth), setup= ''from __main__ import main'', number=N)*1000./N print(template.format(i, depth, n, time)) i += 1 if __name__ == "__main__": main()

La salida de referencia:

i depth n time (ms) 1 4 10 0.182511806488 2 4 50 -- 3 4 100 -- 4 4 150 -- 5 8 10 0.397620201111 6 8 50 1.66054964066 7 8 100 2.90236949921 8 8 150 3.48146915436 9 15 10 0.804011821747 10 15 50 3.7428188324 11 15 100 7.34910964966 12 15 150 10.8230614662 13 16 10 0.804491043091 14 16 50 3.66818904877 15 16 100 7.09567070007 16 16 150 10.404779911 17 17 10 0.865840911865 18 17 50 3.9999294281 19 17 100 7.70257949829 20 17 150 11.3758206367 21 18 10 0.915451049805 22 18 50 4.22935962677 23 18 100 8.22361946106 24 18 150 12.2081303596

10 muestras de tamaño 10 de un árbol con profundidad 10:

[''1111010111'', ''1110111010'', ''1010111010'', ''011110010'', ''0111100001'', ''011101110'', ''01110010'', ''01001111'', ''0001000100'', ''000001010''] [''110'', ''0110101110'', ''0110001100'', ''0011110'', ''0001111011'', ''0001100010'', ''0001100001'', ''0001100000'', ''0000011010'', ''0000001111''] [''11010000'', ''1011111101'', ''1010001101'', ''1001110001'', ''1001100110'', ''10001110'', ''011111110'', ''011001100'', ''0101110000'', ''001110101''] [''11111101'', ''110111'', ''110110111'', ''1101010101'', ''1101001011'', ''1001001100'', ''100100010'', ''0100001010'', ''0100000111'', ''0010010110''] [''111101000'', ''1110111101'', ''1101101'', ''1101000000'', ''1011110001'', ''0111111101'', ''01101011'', ''011010011'', ''01100010'', ''0101100110''] [''1111110001'', ''11000110'', ''1100010100'', ''101010000'', ''1010010001'', ''100011001'', ''100000110'', ''0100001111'', ''001101100'', ''0001101101''] [''111110010'', ''1110100'', ''1101000011'', ''101101'', ''101000101'', ''1000001010'', ''0111100'', ''0101010011'', ''0101000110'', ''000100111''] [''111100111'', ''1110001110'', ''1100111111'', ''1100110010'', ''11000110'', ''1011111111'', ''0111111'', ''0110000100'', ''0100011'', ''0010110111''] [''1101011010'', ''1011111'', ''1011100100'', ''1010000010'', ''10010'', ''1000010100'', ''0111011111'', ''01010101'', ''001101'', ''000101100''] [''111111110'', ''111101001'', ''1110111011'', ''111011011'', ''1001011101'', ''1000010100'', ''0111010101'', ''010100110'', ''0100001101'', ''0010000000'']


Mapeo de ids a cadenas. Puede asignar números a sus vectores 0/1, como mencionó @BrodieG:

# some key objects n_pool = sum(2^(1:max_len)) # total number of indices cuts = cumsum(2^(1:max_len-1)) # new group starts inds_by_g = mapply(seq,cuts,cuts*2) # indices grouped by length # the mapping to strings (one among many possibilities) library(data.table) get_01str <- function(id,max_len){ cuts = cumsum(2^(1:max_len-1)) g = findInterval(id,cuts) gid = id-cuts[g]+1 data.table(g,gid)[,s:= do.call(paste,c(list(sep=""),lapply( seq(g[1]), function(x) (gid-1) %/% 2^(x-1) %% 2 ))) ,by=g]$s }

Encontrar ids para soltar. Vamos a soltar secuencialmente iddesde el grupo de muestreo:

# the mapping from one index to indices of nixed strings get_nixstrs <- function(g,gid,max_len){ cuts = cumsum(2^(1:max_len-1)) gids_child = { x = gid%%2^sequence(g-1) ifelse(x,x,2^sequence(g-1)) } ids_child = gids_child+cuts[sequence(g-1)]-1 ids_parent = if (g==max_len) gid+cuts[g]-1 else { gids_par = vector(mode="list",max_len) gids_par[[g]] = gid for (gg in seq(g,max_len-1)) gids_par[[gg+1]] = c(gids_par[[gg]],gids_par[[gg]]+2^gg) unlist(mapply(`+`,gids_par,cuts-1)) } c(ids_child,ids_parent) }

Los índices están agrupados por g, el número de caracteres nchar(get_01str(id)),. Debido a que los índices están ordenados por g, g=findInterval(id,cuts)es una ruta más rápida.

Un índice en el grupo g, 1 < g < max_lentiene un índice de "niño" de tamaño g-1y dos índices de matrices de tamaño g+1. Para cada nodo hijo, tomamos su nodo hijo hasta que golpeamos g==1; y para cada nodo principal, tomamos su par de nodos principales hasta que pulsamos g==max_len.

La estructura del árbol es la más simple en términos del identificador dentro del grupo gid,. gidmapas a los dos padres, gidy gid+2^g; e invirtiendo este mapeo uno encuentra al niño.

Muestreo

drawem <- function(n,max_len){ cuts = cumsum(2^(1:max_len-1)) inds_by_g = mapply(seq,cuts,cuts*2) oklens = (1:max_len)[ n <= 2^max_len*(1-2^(-(1:max_len)))+1 ] okinds = unlist(inds_by_g[oklens]) mysamp = rep(0,n) for (i in 1:n){ id = if (length(okinds)==1) okinds else sample(okinds,1) g = findInterval(id,cuts) gid = id-cuts[g]+1 nixed = get_nixstrs(g,gid,max_len) # print(id); print(okinds); print(nixed) mysamp[i] = id okinds = setdiff(okinds,nixed) if (!length(okinds)) break } res <- rep("",n) res[seq.int(i)] <- get_01str(mysamp[seq.int(i)],max_len) res }

La oklensparte integra la idea del OP para omitir cadenas que están garantizadas para hacer imposible el muestreo. Sin embargo, incluso haciendo esto, podemos seguir un camino de muestreo que nos deja sin más opciones. Tomando el ejemplo de la OP de max_len=4y n=10, sabemos que hay que bajar 0y 1de la consideración, pero ¿Qué pasa si nuestro primer cuatro empates son 00, 01, 11y 10? Oh bueno, supongo que estamos fuera de suerte. Esta es la razón por la que debería definir las probabilidades de muestreo. (El OP tiene otra idea, para determinar, en cada paso, qué nodos conducirán a un estado imposible, pero eso parece una tarea difícil).

Ilustración

# how the indices line up n_pool = sum(2^(1:max_len)) pdt <- data.table(id=1:n_pool) pdt[,g:=findInterval(id,cuts)] pdt[,gid:=1:.N,by=g] pdt[,s:=get_01str(id,max_len)] # example run set.seed(4); drawem(5,5) # [1] "01100" "1" "0001" "0101" "00101" set.seed(4); drawem(8,4) # [1] "1100" "0" "111" "101" "1101" "100" "" ""

Puntos de referencia (más antiguos que los de la respuesta de @BrodieG)

require(rbenchmark) max_len = 8 n = 8 benchmark( jos_lp = { pool <- unlist(lapply(seq_len(max_len), function(x) do.call(paste0, expand.grid(rep(list(c(''0'', ''1'')), x))))) sample.lp(pool, n)}, bro_string = {pool <- make_pool(max_len);sample01(pool,n)}, fra_num = drawem(n,max_len), replications=5)[1:5] # test replications elapsed relative user.self # 2 bro_string 5 0.05 2.5 0.05 # 3 fra_num 5 0.02 1.0 0.02 # 1 jos_lp 5 1.56 78.0 1.55 n = 12 max_len = 12 benchmark( bro_string={pool <- make_pool(max_len);sample01(pool,n)}, fra_num=drawem(n,max_len), replications=5)[1:5] # test replications elapsed relative user.self # 1 bro_string 5 0.54 6.75 0.51 # 2 fra_num 5 0.08 1.00 0.08

Otras respuestas. Hay otras dos respuestas:

jos_enum = {pool <- unlist(lapply(seq_len(max_len), function(x) do.call(paste0, expand.grid(rep(list(c(''0'', ''1'')), x))))) get.template(pool, n)} bro_num = sample011(max_len,n)

Dejé de lado el método de enumeración de @josilber porque tomó demasiado tiempo; y el método numérico / índice de @ BrodieG porque no funcionó en ese momento, pero ahora sí. Vea la respuesta actualizada de @BrodieG para más evaluaciones comparativas.

Velocidad vs corrección. Si bien las respuestas de @josilber son mucho más lentas (y para el método de enumeración, al parecer mucho más intensivas en memoria), están garantizadas para obtener una muestra de tamaño nen el primer intento. Con el método de cuerdas de @BrodieG o esta respuesta, tendrás que volver a muestrear una y otra vez con la esperanza de dibujar un completo n. Con grandes max_len, eso no debería ser un problema, supongo.

Esta respuesta se escala mejor que bro_stringporque no requiere la construcción del poolfrente.


Un enfoque sería simplemente generar todas las tuplas posibles del tamaño apropiado con un enfoque iterativo:

  1. Construye todas las tuplas de tamaño 1 (todos los elementos en pool)
  2. Tome un producto cruzado con elementos en pool
  3. Remueve cualquier tupla usando el mismo elemento de poolmás de una vez.
  4. Eliminar cualquier duplicado exacto de otra tupla
  5. Retire cualquier tupla con un par que no pueda usarse juntos
  6. Enjuague y repita hasta que obtenga el tamaño de tupla apropiado

Esto es ejecutable para los tamaños dados ( poolde longitud 30, max_len4):

get.template <- function(pool, max_len) { banned <- which(outer(paste0(''^'', pool), pool, Vectorize(grepl)), arr.ind=T) banned <- banned[banned[,1] != banned[,2],] banned <- paste(banned[,1], banned[,2]) vals <- matrix(seq(length(pool))) for (k in 2:max_len) { vals <- cbind(vals[rep(1:nrow(vals), each=length(pool)),], rep(1:length(pool), nrow(vals))) # Can''t sample same value more than once vals <- vals[apply(vals, 1, function(x) length(unique(x)) == length(x)),] # Sort rows to ensure unique only vals <- t(apply(vals, 1, sort)) vals <- unique(vals) # Can''t have banned pair combos <- combn(ncol(vals), 2) for (k in seq(ncol(combos))) { c1 <- combos[1,k] c2 <- combos[2,k] vals <- vals[!paste(vals[,c1], vals[,c2]) %in% banned,] } } return(matrix(pool[vals], nrow=nrow(vals))) } max_len <- 4 pool <- unlist(lapply(seq_len(max_len), function(x) do.call(paste0, expand.grid(rep(list(c(''0'', ''1'')), x))))) system.time(template <- get.template(pool, 4)) # user system elapsed # 4.549 0.050 4.614

Ahora puede muestrear tantas veces como desee de las filas de template(que será muy rápido), y esto será lo mismo que muestrear al azar del espacio definido.


Si no desea generar el conjunto de todas las tuplas posibles y luego muestrearlas de forma aleatoria (lo que, como puede observar, puede ser imposible para tamaños de entrada grandes), otra opción es dibujar una muestra única con programación de enteros. Básicamente, podría asignar cada elemento en poolun valor aleatorio y luego seleccionar la tupla factible con la mayor suma de valores. Esto debería dar a cada tupla la misma probabilidad de ser seleccionado porque todos son del mismo tamaño y sus valores fueron seleccionados al azar. Las restricciones del modelo garantizarían que ninguno de los pares de tuplas no permitidos estén seleccionados y que se seleccione el número correcto de elementos.

Aquí hay una solución con el lpSolvepaquete:

library(lpSolve) sample.lp <- function(pool, max_len) { pool <- sort(pool) pml <- max(nchar(pool)) runs <- c(rev(cumsum(2^(seq(pml-1)))), 0) banned.from <- rep(seq(pool), runs[nchar(pool)]) banned.to <- banned.from + unlist(lapply(runs[nchar(pool)], seq_len)) banned.constr <- matrix(0, nrow=length(banned.from), ncol=length(pool)) banned.constr[cbind(seq(banned.from), banned.from)] <- 1 banned.constr[cbind(seq(banned.to), banned.to)] <- 1 mod <- lp(direction="max", objective.in=runif(length(pool)), const.mat=rbind(banned.constr, rep(1, length(pool))), const.dir=c(rep("<=", length(banned.from)), "=="), const.rhs=c(rep(1, length(banned.from)), max_len), all.bin=TRUE) pool[which(mod$solution == 1)] } set.seed(144) pool <- unlist(lapply(seq_len(4), function(x) do.call(paste0, expand.grid(rep(list(c(''0'', ''1'')), x))))) sample.lp(pool, 4) # [1] "0011" "010" "1000" "1100" sample.lp(pool, 8) # [1] "0000" "0100" "0110" "1001" "1010" "1100" "1101" "1110"

Esto parece escalar a piscinas bastante grandes. Por ejemplo, se tarda un poco más de 2 segundos en obtener una muestra de longitud 20 de un grupo de tamaño 510:

pool <- unlist(lapply(seq_len(8), function(x) do.call(paste0, expand.grid(rep(list(c(''0'', ''1'')), x))))) length(pool) # [1] 510 system.time(sample.lp(pool, 20)) # user system elapsed # 0.232 0.008 0.239

Si realmente necesita resolver los problemas de los problemas realmente grandes, puede pasar de los solucionadores de código no abierto que se envían lpSolvea un solucionador comercial como gurobi o cplex (no es gratis en general, pero es gratuito para uso académico).