Repita las copias de los elementos de la matriz: deencoding de longitud de ejecución en MATLAB

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 cumsumcumsum “, 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 sum 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 XX 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 blog de UNDOCUMENTED MATLAB sobre una preasignación más rápida , 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 desencoding 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

enter image description here

enter image description here

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.

Puntos de referencia

Actualizado para R2015b : repelem ahora más rápido para todos los tamaños de datos.


Funciones probadas:

  1. Función repelem incorporada de repelem que se agregó en R2015a
  2. solución cumsum de cumsum ( rld_cumsum )
  3. La solución cumsum + diff Divakar ( rld_cumsum_diff )
  4. La solución de accumarray de accumarray ( knedlsepp5cumsumccumarray ) de esta publicación
  5. Implementación ingenua basada en bucle ( naive_jit_test.m ) para probar el comstackdor just-in-time

Resultados de test_rld.m en R2015 b :

tiempo repelente

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 que rld_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 que repelem alrededor de 5 000 elementos
  • repelem vuelve más lento que rld_cumsum algún lugar entre 30 000 y 300 000 elementos
  • rld_cumsum tiene aproximadamente el mismo rendimiento que knedlsepp5cumsumccumarray
  • naive_jit_test.m tiene velocidad casi constante y está a la par con rld_cumsum y knedlsepp5cumsumccumarray para tamaños más pequeños, un poco más rápido para tamaños grandes

enter image description here

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 .

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 sum de todas las repeticiones en b ). A continuación, se colocan en el primer elemento y cada elemento subsiguiente representa el lugar en el que se iniciará una nueva secuencia de valores en la salida. La sum 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{:}]; 

Finalmente (a partir de R2015a ) una función incorporada y documentada para hacer esto, repelem . La siguiente syntax, donde el segundo argumento es un vector, es relevante aquí:

W = repelem(V,N) , con el vector V y el vector N , crea un vector W donde el elemento V(i) se repite N(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 progtwig 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:

Parcela actualizada que compara el tiempo de ejecución de repelentes de diferentes métodosSpeedup de repelem sobre cumsum + diff (alrededor de un factor 2)

    Intereting Posts