Generate 1d numpy with pieces of arbitrary length

Question

I need to create a 1D array where repeating sequences of integers are separated by a random number of zeros.

So far I have been using the following code for this:

from random import normalvariate regular_sequence = np.array([1,2,3,4,5], dtype=np.int) n_iter = 10 lag_mean = 10 # mean length of zeros sequence lag_sd = 1 # standard deviation of zeros sequence length # Sequence of lags lengths lag_seq = [int(round(normalvariate(lag_mean, lag_sd))) for x in range(n_iter)] # Generate list of concatenated zeros and regular sequences seq = [np.concatenate((np.zeros(x, dtype=np.int), regular_sequence)) for x in lag_seq] seq = np.concatenate(seq) 

It works, but it looks very slow when I need a lot of long sequences. So how can I optimize it?

Answer 1

You can pre-compute the indexes into which repeating regular_sequence elements should be added, and then set them using regular_sequence in vector form. To pre-compute these indices, you can use np.cumsum to get the beginning of each regular_sequence , and then add a continuous set of integers expanding to the size of regular_sequence to get all the indices you need to update. Thus, the implementation will look something like this:

 # Size of regular_sequence N = regular_sequence.size # Use cumsum to pre-compute start of every occurance of regular_sequence offset_arr = np.cumsum(lag_seq) idx = np.arange(offset_arr.size)*N + offset_arr # Setup output array out = np.zeros(idx.max() + N,dtype=regular_sequence.dtype) # Broadcast the start indices to include entire length of regular_sequence # to get all positions where regular_sequence elements are to be set np.put(out,idx[:,None] + np.arange(N),regular_sequence) 

---

Runtime Tests -

 def original_app(lag_seq, regular_sequence): seq = [np.concatenate((np.zeros(x, dtype=np.int), regular_sequence)) for x in lag_seq] return np.concatenate(seq) def vectorized_app(lag_seq, regular_sequence): N = regular_sequence.size offset_arr = np.cumsum(lag_seq) idx = np.arange(offset_arr.size)*N + offset_arr out = np.zeros(idx.max() + N,dtype=regular_sequence.dtype) np.put(out,idx[:,None] + np.arange(N),regular_sequence) return out In [64]: # Setup inputs ...: regular_sequence = np.array([1,2,3,4,5], dtype=np.int) ...: n_iter = 1000 ...: lag_mean = 10 # mean length of zeros sequence ...: lag_sd = 1 # standard deviation of zeros sequence length ...: ...: # Sequence of lags lengths ...: lag_seq = [int(round(normalvariate(lag_mean, lag_sd))) for x in range(n_iter)] ...: In [65]: out1 = original_app(lag_seq, regular_sequence) In [66]: out2 = vectorized_app(lag_seq, regular_sequence) In [67]: %timeit original_app(lag_seq, regular_sequence) 100 loops, best of 3: 4.28 ms per loop In [68]: %timeit vectorized_app(lag_seq, regular_sequence) 1000 loops, best of 3: 294 µs per loop 

Answer 2

The best approach, I think, would be to use convolution. You can determine the length of the delay, combine this with the length of the sequence, and use it to determine the starting point of each regular sequence. Set these start points to zero, then collapse with your usual sequence to fill in the values.

 import numpy as np regular_sequence = np.array([1,2,3,4,5], dtype=np.int) n_iter = 10000000 lag_mean = 10 # mean length of zeros sequence lag_sd = 1 # standard deviation of zeros sequence length # Sequence of lags lengths lag_lens = np.round(np.random.normal(lag_mean, lag_sd, n_iter)).astype(np.int) lag_lens[1:] += len(regular_sequence) starts_inds = lag_lens.cumsum()-1 # Generate list of convolved ones and regular sequences seq = np.zeros(lag_lens.sum(), dtype=np.int) seq[starts_inds] = 1 seq = np.convolve(seq, regular_sequence) 

This approach takes about 1/20 times in large sequences, even after changing your version to use the numpy random number generator.

Answer 3

This is not a trivial problem, because the data is biased. Performance depends on a long sequence. Let's take an example of a square problem: many, long, regular and zero sequences ( n_iter==n_reg==lag_mean ):

 import numpy as np n_iter = 1000 n_reg = 1000 regular_sequence = np.arange(n_reg, dtype=np.int) lag_mean = n_reg # mean length of zeros sequence lag_sd = lag_mean/10 # standard deviation of zeros sequence length lag_seq=np.int64(np.random.normal(lag_mean,lag_sd,n_iter)) # Sequence of lags lengths 

First your solution:

 def seq_hybrid(): seqs = [np.concatenate((np.zeros(x, dtype=np.int), regular_sequence)) for x in lag_seq] seq = np.concatenate(seqs) return seq 

Then pure numpy one:

 def seq_numpy(): seq=np.zeros(lag_seq.sum()+n_iter*n_reg,dtype=int) cs=np.cumsum(lag_seq+n_reg)-n_reg indexes=np.add.outer(cs,np.arange(n_reg)) seq[indexes]=regular_sequence return seq 

A to solve the cycle:

 def seq_python(): seq=np.empty(lag_seq.sum()+n_iter*n_reg,dtype=int) i=0 for lag in lag_seq: for k in range(lag): seq[i]=0 i+=1 for k in range(n_reg): seq[i]=regular_sequence[k] i+=1 return seq 

And just compilation in time with numba:

 from numba import jit seq_numba=jit(seq_python) 

Tests now:

 In [96]: %timeit seq_hybrid() 10 loops, best of 3: 38.5 ms per loop In [97]: %timeit seq_numpy() 10 loops, best of 3: 34.4 ms per loop In [98]: %timeit seq_python() 1 loops, best of 3: 1.56 s per loop In [99]: %timeit seq_numba() 100 loops, best of 3: 12.9 ms per loop 

In this case, your hybrid solution will be as fast as purely numpy, because performance depends mainly on the inner loop. And yours (zeros and concatenates) is numerical. As expected, the python solution is slower with a traditional 40x ratio. But numpy is not optimal here because it uses the fantastic indexing needed with inconsistent data. In this case, numba can help you: minimal operations are performed at the C level, for a gain of 120x this time compared to the python solution.

For other values ​​of n_iter,n_reg gain compared to python solution:

 n_iter= 1000, n_reg= 1000 : seq_numba 124, seq_hybrid 49, seq_numpy 44. n_iter= 10, n_reg= 100000 : seq_numba 123, seq_hybrid 104, seq_numpy 49. n_iter= 100000, n_reg= 10 : seq_numba 127, seq_hybrid 1, seq_numpy 42. 

Answer 4

I thought the answer posted on this question had a good approach using the binary mask and np.convolve, but the answer was deleted and I don't know why. This is about 2 issues.

 def insert_sequence(lag_seq, regular_sequence): offsets = np.cumsum(lag_seq) start_locs = np.zeros(offsets[-1] + 1, dtype=regular_sequence.dtype) start_locs[offsets] = 1 return np.convolve(start_locs, regular_sequence) lag_seq = np.random.normal(15,1,10) lag_seq = lag_seq.astype(np.uint8) regular_sequence = np.arange(1, 6) seq = insert_sequence(lag_seq, regular_sequence) print(repr(seq))