arrays - Repita las copias de los elementos de la matriz: decodificación de longitud de ejecución en MATLAB
vectorization repeating (5)
Estoy tratando de insertar valores múltiples en una matriz usando una matriz de ''valores'' y una matriz de ''contador''. Por ejemplo, si:
a=[1,3,2,5]
b=[2,2,1,3]
Quiero la salida de alguna función
c=somefunction(a,b)
ser
c=[1,1,3,3,2,5,5,5]
Donde a (1) se repite b (1) número de veces, a (2) recurre b (2) veces, etc ...
¿Hay una función incorporada en MATLAB que hace esto? Me gustaría evitar el uso de un bucle for si es posible. He intentado variaciones de ''repmat ()'' y ''kron ()'' en vano.
Esto es básicamente Run-length encoding
.
Planteamiento del problema
Tenemos una variedad de valores, vals
y runlengths, runlens
:
vals = [1,3,2,5]
runlens = [2,2,1,3]
Necesitamos repetir cada elemento en vals
cada elemento correspondiente en los runlens
. Por lo tanto, el resultado final sería:
output = [1,1,3,3,2,5,5,5]
Enfoque prospectivo
Una de las herramientas más rápidas con MATLAB es cumsum
y es muy útil cuando se trata de vectorizar problemas que funcionan en patrones irregulares. En el problema indicado, la irregularidad viene con los diferentes elementos en los runlens
.
Ahora, para explotar cumsum
, tenemos que hacer dos cosas aquí: Inicializar una matriz de zeros
y colocar valores "apropiados" en posiciones "clave" sobre la matriz de ceros, de modo que después de cumsum
" cumsum
", terminemos con una conjunto final de vals
repetidos de tiempos runlens
.
Pasos: enumeremos los pasos mencionados anteriormente para dar una perspectiva más fácil al enfoque prospectivo:
1) Initialize zeros array: ¿Cuál debe ser la longitud? Dado que estamos repitiendo tiempos de runlens
, la longitud de la matriz de ceros debe ser la suma de todos los runlens
.
2) Buscar posiciones / índices clave: ahora estas posiciones clave son lugares a lo largo de la matriz de ceros donde cada elemento de vals
comienza a repetirse. Por lo tanto, para runlens = [2,2,1,3]
, las posiciones clave asignadas a la matriz de ceros serían:
[X 0 X 0 X X 0 0], where X''s are those key positions.
3) Encuentre valores apropiados: el último clavo que se martillará antes de usar cumsum
sería poner valores "apropiados" en esas posiciones clave. Ahora, ya que estaríamos haciendo cumsum
poco después, si lo piensas bien, necesitarías una versión differentiated
de values
con diff
, para que cumsum
sobre esos devuelva nuestros values
. Dado que estos valores diferenciados se colocarían en una matriz de ceros en lugares separados por las distancias de los runlens
, después de usar cumsum
tendríamos cada elemento de runlens
veces runlens
repetidos como salida final.
Código de solución
Aquí está la implementación que une todos los pasos mencionados anteriormente:
%// Calculate cumsumed values of runLengths.
%// We would need this to initialize zeros array and find key positions later on.
clens = cumsum(runlens)
%// Initalize zeros array
array = zeros(1,(clens(end)))
%// Find key positions/indices
key_pos = [1 clens(1:end-1)+1]
%// Find appropriate values
app_vals = diff([0 vals])
%// Map app_values at key_pos on array
array(pos) = app_vals
%// cumsum array for final output
output = cumsum(array)
Hack de preasignación
Como se puede ver, el código mencionado anteriormente usa preasignación con ceros. Ahora, de acuerdo con este UNDOCUMENTED MATLAB blog on faster pre-allocation
, se puede lograr una preasignación mucho más rápida con -
`array(clens(end)) = 0` instead of `array = zeros(1,(clens(end)))`
Embalaje: código de función
Para concluir todo, tendríamos un código de función compacto para lograr esta descodificación de longitud de ejecución como tal -
function out = rle_cumsum_diff(vals,runlens)
clens = cumsum(runlens);
idx(clens(end))=0;
idx([1 clens(1:end-1)+1]) = diff([0 vals]);
out = cumsum(idx);
return;
Benchmarking
Código de Benchmarking
A continuación se enumera el código de evaluación comparativa para comparar tiempos de ejecución y aceleraciones para el enfoque cumsum+diff
en esta publicación sobre el otro enfoque basado en MATLAB 2014B
cumsum-only
en MATLAB 2014B
-
datasizes = [reshape(linspace(10,70,4).''*10.^(0:4),1,[]) 10^6 2*10^6]; %//''
fcns = {''rld_cumsum'',''rld_cumsum_diff''}; %// approaches to be benchmarked
for k1 = 1:numel(datasizes)
n = datasizes(k1); %// Create random inputs
vals = randi(200,1,n);
runs = [5000 randi(200,1,n-1)]; %// 5000 acts as an aberration
for k2 = 1:numel(fcns) %// Time approaches
tsec(k2,k1) = timeit(@() feval(fcns{k2}, vals,runs), 1);
end
end
figure, %// Plot runtimes
loglog(datasizes,tsec(1,:),''-bo''), hold on
loglog(datasizes,tsec(2,:),''-k+'')
set(gca,''xgrid'',''on''),set(gca,''ygrid'',''on''),
xlabel(''Datasize ->''), ylabel(''Runtimes (s)'')
legend(upper(strrep(fcns,''_'','' ''))),title(''Runtime Plot'')
figure, %// Plot speedups
semilogx(datasizes,tsec(1,:)./tsec(2,:),''-rx'')
set(gca,''ygrid'',''on''), xlabel(''Datasize ->'')
legend(''Speedup(x) with cumsum+diff over cumsum-only''),title(''Speedup Plot'')
Código de función asociada para rld_cumsum.m
:
function out = rld_cumsum(vals,runlens)
index = zeros(1,sum(runlens));
index([1 cumsum(runlens(1:end-1))+1]) = 1;
out = vals(cumsum(index));
return;
Runtime y Speedup Parcelas
Conclusiones
El enfoque propuesto parece estar dándonos una aceleración notable sobre el cumsum-only
, ¡que es aproximadamente 3x
!
¿Por qué es este nuevo cumsum+diff
mejor que el enfoque anterior de cumsum-only
?
Bueno, la esencia de la razón se encuentra en el paso final del cumsum-only
que necesita para mapear los valores "cumsumed" en vals
. En el nuevo cumsum+diff
, estamos haciendo diff(vals)
para el cual MATLAB procesa solo n
elementos (donde n es el número de longitudes de ejecución) en comparación con el mapeo de sum(runLengths)
número de elementos para el cumsum-only
enfoque y este número debe ser muchas veces más que n
por lo tanto, la notable aceleración con este nuevo enfoque.
Finalmente (a partir de R2015a ) una función incorporada y documentada para hacer esto, repelem
. La siguiente sintaxis, donde el segundo argumento es un vector, es relevante aquí:
W = repelem(V,N)
, con el vectorV
y el vectorN
, crea un vectorW
donde el elementoV(i)
se repiteN(i)
veces.
O dicho de otra manera, "Cada elemento de N
especifica el número de veces que se debe repetir el elemento correspondiente de V
".
Ejemplo:
>> a=[1,3,2,5]
a =
1 3 2 5
>> b=[2,2,1,3]
b =
2 2 1 3
>> repelem(a,b)
ans =
1 1 3 3 2 5 5 5
Los problemas de rendimiento en el repelem
integrado de repelem
se han corregido a partir del R2015b. He ejecutado el programa test_rld.m
desde la publicación de chappjc en R2015b, y repelem
ahora es más rápido que otros algoritmos en un factor 2 aproximadamente:
No conozco ninguna función incorporada, pero aquí hay una solución:
index = zeros(1,sum(b));
index([1 cumsum(b(1:end-1))+1]) = 1;
c = a(cumsum(index));
Explicación:
Primero se crea un vector de ceros de la misma longitud que la matriz de salida (es decir, la suma de todas las repeticiones en b
). Luego, se colocan los elementos en el primer elemento y cada elemento subsiguiente representa el lugar donde comenzará el inicio de una nueva secuencia de valores en la salida. La suma acumulativa del index
vectorial se puede usar para indexar en a
, replicando cada valor el número deseado de veces.
En aras de la claridad, así es como se ven los diversos vectores para los valores de a
y b
dados en la pregunta:
index = [1 0 1 0 1 1 0 0]
cumsum(index) = [1 1 2 2 3 4 4 4]
c = [1 1 3 3 2 5 5 5]
EDITAR: en aras de la exhaustividad, hay otra alternativa que usa ARRAYFUN , pero parece que lleva entre 20 y 100 veces más tiempo ejecutar que la solución anterior con vectores de hasta 10.000 elementos de longitud:
c = arrayfun(@(x,y) x.*ones(1,y),a,b,''UniformOutput'',false);
c = [c{:}];
Puntos de referencia
Actualizado para R2015b : repelem
ahora más rápido para todos los tamaños de datos.
Funciones probadas:
- Función
repelem
incorporada derepelem
que se agregó en R2015a - solución cumsum de
cumsum
(rld_cumsum
) - La solución
cumsum
+diff
Divakar (rld_cumsum_diff
) - La solución de accumarray de
accumarray
(knedlsepp5cumsumaccumarray
) de esta publicación - Implementación ingenua basada en bucle (
naive_jit_test.m
) para probar el compilador just-in-time
Resultados de test_rld.m
en R2015 b :
Antigua gráfica de tiempos usando R2015 a aquí .
Hallazgos :
-
repelem
es siempre el más rápido por aproximadamente un factor de 2. -
rld_cumsum_diff
es consistentemente más rápido querld_cumsum
. -
repelem
es más rápido para tamaños de datos pequeños (menos de 300-500 elementos) -
rld_cumsum_diff
vuelve significativamente más rápido querepelem
alrededor de 5 000 elementos -
repelem
vuelve más lento querld_cumsum
algún lugar entre 30 000 y 300 000 elementos -
rld_cumsum
tiene aproximadamente el mismo rendimiento queknedlsepp5cumsumaccumarray
-
naive_jit_test.m
tiene velocidad casi constante y está a la par conrld_cumsum
yknedlsepp5cumsumaccumarray
para tamaños más pequeños, un poco más rápido para tamaños grandes
Antigua gráfica de tasas usando R2015 a aquí .
Conclusión
Use repelem
debajo de unos 5 000 elementos y la solución . cumsum
+ diff
anterior