python - structures - La subsecuencia más espaciada equitativamente
the algorithms (10)
Tengo un millón de enteros ordenados y me gustaría encontrar la subsecuencia más larga donde la diferencia entre pares consecutivos es la misma. Por ejemplo
1, 4, 5, 7, 8, 12
tiene una subsecuencia
4, 8, 12
Mi método ingenuo es codicioso y solo comprueba hasta dónde puedes extender una subsecuencia desde cada punto. Esto toma O(n²)
tiempo por punto que parece.
¿Hay una forma más rápida de resolver este problema?
Actualizar. Voy a probar el código que figura en las respuestas lo antes posible (gracias). Sin embargo, ya está claro que el uso de la memoria n ^ 2 no funcionará. Hasta el momento no hay un código que finalice con la entrada como [random.randint(0,100000) for r in xrange(200000)]
.
Tiempos Probé con los siguientes datos de entrada en mi sistema de 32 bits.
a= [random.randint(0,10000) for r in xrange(20000)]
a.sort()
- El método de programación dinámica de ZelluX usa 1.6G de RAM y toma 2 minutos y 14 segundos. ¡Con pypy solo lleva 9 segundos! Sin embargo, falla con un error de memoria en entradas grandes.
- El método O (nd) time de Armin tomó 9 segundos con pypy pero solo 20MB de RAM. Por supuesto, esto sería mucho peor si el alcance fuera mucho mayor. El bajo uso de memoria significaba que también podía probarlo con a = [random.randint (0,100000) para r en xrange (200000)] pero no terminó en los pocos minutos que se lo di con pypy.
Para poder probar el método de Kluev, reran con
a= [random.randint(0,40000) for r in xrange(28000)]
a = list(set(a))
a.sort()
para hacer una lista de longitud aproximadamente 20000
. Todos los tiempos con pypy
- ZelluX, 9 segundos
- Kluev, 20 segundos
- Armin, 52 segundos
Parece que si el método ZelluX pudiera hacerse espacio lineal, sería el claro ganador.
ACTUALIZACIÓN: He encontrado un documento sobre este problema, puedes descargarlo here .
Aquí hay una solución basada en programación dinámica. Requiere O (n ^ 2) complejidad de tiempo y O (n ^ 2) complejidad de espacio, y no usa hash.
Suponemos que todos los números se guardan en el conjunto a
en orden ascendente, y n
guarda su longitud. La matriz 2D l[i][j]
define la longitud de la subsecuencia más espaciada igualmente equidistante que termina con a[i]
y a[j]
, y l[j][k]
= l[i][j]
+ 1 si a[j]
- a[i]
= a[k]
- a[j]
(i <j <k).
lmax = 2
l = [[2 for i in xrange(n)] for j in xrange(n)]
for mid in xrange(n - 1):
prev = mid - 1
succ = mid + 1
while (prev >= 0 and succ < n):
if a[prev] + a[succ] < a[mid] * 2:
succ += 1
elif a[prev] + a[succ] > a[mid] * 2:
prev -= 1
else:
l[mid][succ] = l[prev][mid] + 1
lmax = max(lmax, l[mid][succ])
prev -= 1
succ += 1
print lmax
Aquí hay otra respuesta, que funciona en el tiempo O(n^2)
y sin requisitos de memoria notables más allá de convertir la lista en un conjunto.
La idea es bastante ingenua: al igual que el póster original, es codicioso y solo comprueba hasta dónde puedes extender una subsecuencia de cada par de puntos, sin embargo, verificando primero que estamos al comienzo de una subsecuencia. En otras palabras, desde los puntos a
y b
compruebas qué tan lejos puedes extender a b + (ba)
, b + 2*(ba)
, ... pero solo si a - (ba)
no está ya en el conjunto de todos puntos. Si es así, entonces ya viste la misma subsecuencia.
El truco consiste en convencernos de que esta simple optimización es suficiente para reducir la complejidad a O(n^2)
del O(n^3)
original O(n^3)
. Eso queda como un ejercicio para el lector :-) El tiempo es competitivo con otras soluciones O(n^2)
aquí.
A = [1, 4, 5, 7, 8, 12] # in sorted order
Aset = set(A)
lmax = 2
for j, b in enumerate(A):
for i in range(j):
a = A[i]
step = b - a
if b + step in Aset and a - step not in Aset:
c = b + step
count = 3
while c + step in Aset:
c += step
count += 1
#print "found %d items in %d .. %d" % (count, a, c)
if count > lmax:
lmax = count
print lmax
Este es mi 2 centavos.
Si tiene una lista llamada entrada:
input = [1, 4, 5, 7, 8, 12]
Puede construir una estructura de datos que para cada uno de estos puntos (excluyendo el primero) le dirá qué tan lejos está ese punto de cualquiera de sus predecesores:
[1, 4, 5, 7, 8, 12]
x 3 4 6 7 11 # distance from point i to point 0
x x 1 3 4 8 # distance from point i to point 1
x x x 2 3 7 # distance from point i to point 2
x x x x 1 5 # distance from point i to point 3
x x x x x 4 # distance from point i to point 4
Ahora que tiene las columnas, puede considerar el i-th
elemento de entrada (que es la input[i]
) y cada número n
en su columna.
Los números que pertenecen a una serie de números equidistantes que incluyen la input[i]
, son aquellos que tienen n * j
en la posición i-th
de su columna, donde j
es el número de coincidencias encontradas al mover las columnas de izquierda a derecha , más el k-th
predecesor de la input[i]
, donde k
es el índice de n
en la columna de input[i]
.
Ejemplo: si consideramos i = 1
, input[i] = 4
, n = 3
, entonces, podemos identificar una secuencia que comprenda 4
( input[i]
), 7
(porque tiene un 3
en la posición 1
de su columna) y 1
, porque k
es 0, entonces tomamos el primer predecesor de i
.
Posible implementación (lo siento si el código no usa la misma notación que la explicación):
def build_columns(l):
columns = {}
for x in l[1:]:
col = []
for y in l[:l.index(x)]:
col.append(x - y)
columns[x] = col
return columns
def algo(input, columns):
seqs = []
for index1, number in enumerate(input[1:]):
index1 += 1 #first item was sliced
for index2, distance in enumerate(columns[number]):
seq = []
seq.append(input[index2]) # k-th pred
seq.append(number)
matches = 1
for successor in input[index1 + 1 :]:
column = columns[successor]
if column[index1] == distance * matches:
matches += 1
seq.append(successor)
if (len(seq) > 2):
seqs.append(seq)
return seqs
El más largo:
print max(sequences, key=len)
Este es un caso particular para el problema más genérico que se describe aquí: Descubrir patrones largos donde K = 1 y es fijo. Se demuestra allí que se puede resolver en O (N ^ 2). Ejecutando mi implementación del algoritmo C propuesto, lleva 3 segundos encontrar la solución para N = 20000 y M = 28000 en mi máquina de 32 bits.
Método codicioso
1. Solo se genera una secuencia de decisión.
2. Se generan muchas decisiones. Programación dinámica 1. No garantiza dar una solución óptima siempre.
2. Definitivamente da una solución óptima.
Podemos tener una solución O(n*m)
en el tiempo con muy poca memoria, adaptando la tuya. Aquí n
es el número de elementos en la secuencia de entrada de números dada, y m
es el rango, es decir, el número más alto menos el más bajo.
Llame a A la secuencia de todos los números de entrada (y use un set()
precalculado set()
para responder en tiempo constante a la pregunta "¿este número está en A?"). Llame al paso de la subsecuencia que estamos buscando (la diferencia entre dos números de esta subsecuencia). Para cada valor posible de d, realice la siguiente exploración lineal sobre todos los números de entrada: para cada número n de A en orden creciente, si el número aún no se ha visto, mire hacia adelante en A para la longitud de la secuencia que comienza en n con paso d. A continuación, marque todos los elementos en esa secuencia como ya se ha visto, de modo que evitemos buscarlos nuevamente, por el mismo d. Debido a esto, la complejidad es solo O(n)
para cada valor de d.
A = [1, 4, 5, 7, 8, 12] # in sorted order
Aset = set(A)
for d in range(1, 12):
already_seen = set()
for a in A:
if a not in already_seen:
b = a
count = 1
while b + d in Aset:
b += d
count += 1
already_seen.add(b)
print "found %d items in %d .. %d" % (count, a, b)
# collect here the largest ''count''
Actualizaciones:
Esta solución puede ser lo suficientemente buena si solo te interesan los valores de d que son relativamente pequeños; por ejemplo, si obtener el mejor resultado para
d <= 1000
sería lo suficientemente bueno. Entonces la complejidad baja aO(n*1000)
. Esto hace que el algoritmo sea aproximativo, pero realmente ejecutable paran=1000000
. (Medido a 400-500 segundos con CPython, 80-90 segundos con PyPy, con un subconjunto aleatorio de números entre 0 y 10''000''000).Si aún desea buscar todo el rango, y si el caso común es que existen secuencias largas, una mejora notable es detenerse tan pronto como d sea demasiado grande para encontrar una secuencia aún más larga.
Recorre la matriz, manteniendo un registro de los resultados óptimos y una tabla con
(1) índice - la diferencia del elemento en la secuencia,
(2) recuento: cantidad de elementos en la secuencia hasta el momento, y
(3) el último elemento registrado.
Para cada elemento de matriz, observe la diferencia de cada elemento de matriz anterior; si ese elemento está último en una secuencia indexada en la tabla, ajuste esa secuencia en la tabla y actualice la mejor secuencia, si corresponde, de lo contrario comience una nueva secuencia, a menos que la máxima actual sea mayor que la longitud de la posible secuencia.
Escaneando hacia atrás podemos detener nuestro escaneo cuando d es mayor que el medio del rango de la matriz; o cuando el máximo actual es mayor que la longitud de la secuencia posible, para d mayor que la diferencia indexada más grande. Las secuencias donde s[j]
es mayor que el último elemento de la secuencia se eliminan.
Convertí mi código de JavaScript a Python (mi primer código python):
import random
import timeit
import sys
#s = [1,4,5,7,8,12]
#s = [2, 6, 7, 10, 13, 14, 17, 18, 21, 22, 23, 25, 28, 32, 39, 40, 41, 44, 45, 46, 49, 50, 51, 52, 53, 63, 66, 67, 68, 69, 71, 72, 74, 75, 76, 79, 80, 82, 86, 95, 97, 101, 110, 111, 112, 114, 115, 120, 124, 125, 129, 131, 132, 136, 137, 138, 139, 140, 144, 145, 147, 151, 153, 157, 159, 161, 163, 165, 169, 172, 173, 175, 178, 179, 182, 185, 186, 188, 195]
#s = [0, 6, 7, 10, 11, 12, 16, 18, 19]
m = [random.randint(1,40000) for r in xrange(20000)]
s = list(set(m))
s.sort()
lenS = len(s)
halfRange = (s[lenS-1] - s[0]) // 2
while s[lenS-1] - s[lenS-2] > halfRange:
s.pop()
lenS -= 1
halfRange = (s[lenS-1] - s[0]) // 2
while s[1] - s[0] > halfRange:
s.pop(0)
lenS -=1
halfRange = (s[lenS-1] - s[0]) // 2
n = lenS
largest = (s[n-1] - s[0]) // 2
#largest = 1000 #set the maximum size of d searched
maxS = s[n-1]
maxD = 0
maxSeq = 0
hCount = [None]*(largest + 1)
hLast = [None]*(largest + 1)
best = {}
start = timeit.default_timer()
for i in range(1,n):
sys.stdout.write(repr(i)+"/r")
for j in range(i-1,-1,-1):
d = s[i] - s[j]
numLeft = n - i
if d != 0:
maxPossible = (maxS - s[i]) // d + 2
else:
maxPossible = numLeft + 2
ok = numLeft + 2 > maxSeq and maxPossible > maxSeq
if d > largest or (d > maxD and not ok):
break
if hLast[d] != None:
found = False
for k in range (len(hLast[d])-1,-1,-1):
tmpLast = hLast[d][k]
if tmpLast == j:
found = True
hLast[d][k] = i
hCount[d][k] += 1
tmpCount = hCount[d][k]
if tmpCount > maxSeq:
maxSeq = tmpCount
best = {''len'': tmpCount, ''d'': d, ''last'': i}
elif s[tmpLast] < s[j]:
del hLast[d][k]
del hCount[d][k]
if not found and ok:
hLast[d].append(i)
hCount[d].append(2)
elif ok:
if d > maxD:
maxD = d
hLast[d] = [i]
hCount[d] = [2]
end = timeit.default_timer()
seconds = (end - start)
#print (hCount)
#print (hLast)
print(best)
print(seconds)
Su solución es O(N^3)
ahora (usted dijo O(N^2) per index
). Aquí está O(N^2)
de tiempo y O(N^2)
de solución de memoria.
Idea
Si conocemos la subsecuencia que pasa por los índices i[0]
, i[1]
, i[2]
, i[3]
no deberíamos intentar la subsecuencia que comienza con i[1]
i[2]
o i[2]
y i[3]
Tenga en cuenta que edité ese código para hacerlo un poco más fácil usando eso, pero no funcionará para elementos iguales. Puede verificar fácilmente el número máximo de elementos iguales en O(N)
Pseudocódigo
Solo busco la longitud máxima, pero eso no cambia nada
whereInA = {}
for i in range(n):
whereInA[a[i]] = i; // It doesn''t matter which of same elements it points to
boolean usedPairs[n][n];
for i in range(n):
for j in range(i + 1, n):
if usedPair[i][j]:
continue; // do not do anything. It was in one of prev sequences.
usedPair[i][j] = true;
//here quite stupid solution:
diff = a[j] - a[i];
if diff == 0:
continue; // we can''t work with that
lastIndex = j
currentLen = 2
while whereInA contains index a[lastIndex] + diff :
nextIndex = whereInA[a[lastIndex] + diff]
usedPair[lastIndex][nextIndex] = true
++currentLen
lastIndex = nextIndex
// you may store all indicies here
maxLen = max(maxLen, currentLen)
Pensamientos sobre el uso de la memoria
O(n^2)
es muy lento para 1000000 elementos. Pero si va a ejecutar este código en dicha cantidad de elementos, el mayor problema será el uso de la memoria.
¿Qué se puede hacer para reducirlo?
- Cambie las matrices booleanas a bitfields para almacenar más booleanos por bit.
- Haga que cada siguiente matriz booleana sea más corta porque solo usamos
usedPairs[i][j]
sii < j
Pocas heurísticas:
- Almacene solo pares de indices usados. (Conflictos con la primera idea)
- Elimine los pares usados que nunca usarán más (que son para tal
i
,j
que ya fueron elegidos en el ciclo)
Actualización: el primer algoritmo descrito aquí está obsoleto por la segunda respuesta de Armin Rigo , que es mucho más simple y más eficiente. Pero ambos métodos tienen una desventaja. Necesitan muchas horas para encontrar el resultado de un millón de enteros. Así que probé dos variantes más (ver la segunda mitad de esta respuesta) donde se supone que el rango de los enteros de entrada es limitado. Tal limitación permite algoritmos mucho más rápidos. También traté de optimizar el código de Armin Rigo. Vea mis resultados de evaluación comparativa al final.
Aquí hay una idea del algoritmo que usa la memoria O (N). La complejidad del tiempo es O (N 2 log N), pero puede reducirse a O (N 2 ).
Algoritmo usa las siguientes estructuras de datos:
-
prev
: matriz de índices que apunta al elemento anterior de la subsecuencia (posiblemente incompleta). -
hash
: hashmap con key = diferencia entre pares consecutivos en subsequence y value = otros dos hashmaps. Para estos otros hashmaps: clave = índice de inicio / finalización de la subsecuencia, valor = par de (longitud de subsecuencia, índice de finalización / inicio de la subsecuencia). -
pq
: cola de prioridad para todos los posibles valores de "diferencia" para subsecuencias almacenadas enprev
yhash
.
Algoritmo:
- Inicializar
prev
con los índicesi-1
. Actualicehash
ypq
para registrar todas las subsecuencias (incompletas) encontradas en este paso y sus "diferencias". - Obtenga (y elimine) la "diferencia" más pequeña de
pq
. Obtenga el registro correspondiente delhash
y escanee uno de los mapas hash de segundo nivel. En este momento todas las subsecuencias con "diferencia" dada están completas. Si el mapa hash de segundo nivel contiene una longitud de subsecuencia mejor que la encontrada hasta ahora, actualice el mejor resultado. - En la matriz
prev
: para cada elemento de cualquier secuencia encontrada en el paso 2, disminuya el índice y actualice elhash
y posiblementepq
. Al actualizarhash
, podríamos realizar una de las siguientes operaciones: agregar una subsecuencia nueva de longitud 1, o hacer crecer alguna subsecuencia existente por 1, o fusionar dos subsecuencias existentes. - Elimine el registro del mapa de hash que se encuentra en el paso # 2.
- Continúe desde el paso # 2 mientras
pq
no está vacío.
Este algoritmo actualiza O (N) elementos de prev
O (N) veces cada uno. Y cada una de estas actualizaciones puede requerir agregar una nueva "diferencia" a pq
. Todo esto significa complejidad temporal de O (N 2 log N) si usamos la implementación de pila simple para pq
. Para disminuirlo a O (N 2 ) podríamos usar implementaciones de cola de prioridad más avanzada. Algunas de las posibilidades se enumeran en esta página: Colas de prioridad .
Ver el código Python correspondiente en Ideone . Este código no permite elementos duplicados en la lista. Es posible solucionar esto, pero sería una buena optimización para eliminar duplicados (y para encontrar la subsecuencia más larga más allá de los duplicados por separado).
Y el mismo código después de una pequeña optimización . Aquí la búsqueda finaliza tan pronto como la longitud de la subsecuencia multiplicada por una subsecuencia posible "diferencia" excede el rango de la lista fuente.
El código de Armin Rigo es simple y bastante eficiente. Pero en algunos casos hace algunos cálculos adicionales que pueden evitarse. La búsqueda puede finalizar tan pronto como la longitud de la subsecuencia multiplicada por una posible "diferencia" de subsecuencias exceda el rango de la lista fuente:
def findLESS(A):
Aset = set(A)
lmax = 2
d = 1
minStep = 0
while (lmax - 1) * minStep <= A[-1] - A[0]:
minStep = A[-1] - A[0] + 1
for j, b in enumerate(A):
if j+d < len(A):
a = A[j+d]
step = a - b
minStep = min(minStep, step)
if a + step in Aset and b - step not in Aset:
c = a + step
count = 3
while c + step in Aset:
c += step
count += 1
if count > lmax:
lmax = count
d += 1
return lmax
print(findLESS([1, 4, 5, 7, 8, 12]))
Si el rango de enteros en los datos de origen (M) es pequeño, es posible un algoritmo simple con el tiempo O (M 2 ) y el espacio O (M):
def findLESS(src):
r = [False for i in range(src[-1]+1)]
for x in src:
r[x] = True
d = 1
best = 1
while best * d < len(r):
for s in range(d):
l = 0
for i in range(s, len(r), d):
if r[i]:
l += 1
best = max(best, l)
else:
l = 0
d += 1
return best
print(findLESS([1, 4, 5, 7, 8, 12]))
Es similar al primer método de Armin Rigo, pero no utiliza ninguna estructura de datos dinámica. Supongo que los datos fuente no tienen duplicados. Y (para mantener el código simple), también supongo que el valor mínimo de entrada es no negativo y cercano a cero.
Algoritmo previo puede mejorarse si en lugar de la matriz de booleanos usamos una estructura de datos bitset y operaciones bit a bit para procesar datos en paralelo. El código que se muestra a continuación implementa el conjunto de bits como un número entero integrado de Python. Tiene las mismas suposiciones: no hay duplicados, el valor mínimo de entrada no es negativo y está cerca de cero. La complejidad del tiempo es O (M 2 * log L) donde L es la longitud de la subsecuencia óptima, la complejidad del espacio es O (M):
def findLESS(src):
r = 0
for x in src:
r |= 1 << x
d = 1
best = 1
while best * d < src[-1] + 1:
c = best
rr = r
while c & (c-1):
cc = c & -c
rr &= rr >> (cc * d)
c &= c-1
while c != 1:
c = c >> 1
rr &= rr >> (c * d)
rr &= rr >> d
while rr:
rr &= rr >> d
best += 1
d += 1
return best
Puntos de referencia:
Los datos de entrada (alrededor de 100000 enteros) se generan de esta manera:
random.seed(42)
s = sorted(list(set([random.randint(0,200000) for r in xrange(140000)])))
Y para los algoritmos más rápidos también utilicé los siguientes datos (alrededor de 1000000 enteros):
s = sorted(list(set([random.randint(0,2000000) for r in xrange(1400000)])))
Todos los resultados muestran el tiempo en segundos:
Size: 100000 1000000
Second answer by Armin Rigo: 634 ?
By Armin Rigo, optimized: 64 >5000
O(M^2) algorithm: 53 2940
O(M^2*L) algorithm: 7 711
Algoritmo
- Bucle principal que atraviesa la lista
- Si el número se encuentra en la lista de precalculación, entonces pertenece a todas las secuencias que están en esa lista, recalcula todas las secuencias con conteo + 1
- Eliminar todo precalculado para el elemento actual
- Recalcule nuevas secuencias donde el primer elemento está del rango de 0 a actual, y segundo es el elemento actual de recorrido (en realidad, no de 0 a actual, podemos usar el hecho de que el nuevo elemento no debe ser más que max (a) y nuevo la lista debe tener la posibilidad de ser más larga que la que ya se encontró)
Entonces, para la lista [1, 2, 4, 5, 7]
salida sería (es un poco desordenado, pruebe el código usted mismo y vea)
- índice 0 , elemento 1 :
- si
1
en precalc? No, no hagas nada - Hacer nada
- si
- índice 1 , elemento 2 :
- si
2
en precalc? No, no hagas nada - comprobar si 3 =
1
+ (2
-1
) * 2 en nuestro conjunto? No, no hagas nada
- si
- índice 2 , elemento 4 :
- si
4
en precalc? No, no hagas nada- comprobar si 6 =
2
+ (4
-2
) * 2 en nuestro conjunto? No - comprobar si 7 =
1
+ (4
-1
) * 2 en nuestro conjunto? Sí, agregue un nuevo elemento{7: {3: {''count'': 2, ''start'': 1}}}
7 - elemento de la lista, 3 es paso.
- comprobar si 6 =
- si
- índice 3 , elemento
5
:- si
5
en precalc? No, no hagas nada- no marque
4
porque 6 = 4 + (5
-4
) * 2 es menos que el elemento calculado 7 - comprobar si 8 =
2
+ (5
-2
) * 2 en nuestro conjunto? No - marque 10 =
2
+ (5
-1
) * 2 - más que max (a) == 7
- no marque
- si
- índice 4 , elemento
7
:- si 7 en precalc? Sí, ponlo en el resultado
- no marque
5
porque 9 = 5 + (7
-5
) * 2 es más que max (a) == 7
- no marque
- si 7 en precalc? Sí, ponlo en el resultado
result = (3, {''count'': 3, ''start'': 1}) # step 3, count 3, start 1, conviértalo en secuencia
Complejidad
No debe ser más que O (N ^ 2), y creo que es menor debido a la terminación anticipada de la búsqueda de nuevas secuencias. Trataré de proporcionar un análisis detallado más adelante.
Código
def add_precalc(precalc, start, step, count, res, N):
if step == 0: return True
if start + step * res[1]["count"] > N: return False
x = start + step * count
if x > N or x < 0: return False
if precalc[x] is None: return True
if step not in precalc[x]:
precalc[x][step] = {"start":start, "count":count}
return True
def work(a):
precalc = [None] * (max(a) + 1)
for x in a: precalc[x] = {}
N, m = max(a), 0
ind = {x:i for i, x in enumerate(a)}
res = (0, {"start":0, "count":0})
for i, x in enumerate(a):
for el in precalc[x].iteritems():
el[1]["count"] += 1
if el[1]["count"] > res[1]["count"]: res = el
add_precalc(precalc, el[1]["start"], el[0], el[1]["count"], res, N)
t = el[1]["start"] + el[0] * el[1]["count"]
if t in ind and ind[t] > m:
m = ind[t]
precalc[x] = None
for y in a[i - m - 1::-1]:
if not add_precalc(precalc, y, x - y, 2, res, N): break
return [x * res[0] + res[1]["start"] for x in range(res[1]["count"])]