title: How can I speed up transition matrix creation in Numpy?
body: <p>The following is the most basic way I know of to count transitions in a markov chain and use it to populate a transition matrix:</p>
  
  <pre><code>def increment_counts_in_matrix_from_chain(markov_chain, transition_counts_matrix):
      for i in xrange(1, len(markov_chain)):
          old_state = markov_chain[i - 1]
          new_state = markov_chain[i]
          transition_counts_matrix[old_state, new_state] += 1
  </code></pre>
  
  <p>I've tried speeding it up in 3 different ways:</p>
  
  <p>1) Using a sparse matrix one-liner based on this Matlab code:</p>
  
  <pre><code>transition_matrix = full(sparse(markov_chain(1:end-1), markov_chain(2:end), 1))
  </code></pre>
  
  <p>Which in Numpy/SciPy, looks like this:</p>
  
  <pre><code>def get_sparse_counts_matrix(markov_chain, number_of_states):
      return coo_matrix(([1]*(len(markov_chain) - 1), (markov_chain[0:-1], markov_chain[1:])), shape=(number_of_states, number_of_states)) 
  </code></pre>
  
  <p>And I've tried a couple more Python tweaks, like using zip():</p>
  
  <pre><code>for old_state, new_state in zip(markov_chain[0:-1], markov_chain[1:]):
      transition_counts_matrix[old_state, new_state] += 1 
  </code></pre>
  
  <p>And Queues:</p>
  
  <pre><code>old_and_new_states_holder = Queue(maxsize=2)
  old_and_new_states_holder.put(markov_chain[0])
  for new_state in markov_chain[1:]:
      old_and_new_states_holder.put(new_state)
      old_state = old_and_new_states_holder.get()
      transition_counts_matrix[old_state, new_state] += 1
  </code></pre>
  
  <p>But none of these 3 methods sped things up.  In fact, everything but the zip() solution was at least 10X slower than my original solution.  </p>
  
  <p>Are there any other solutions worth looking into?</p>
  
  <p><br>
  <br>
  <strong>Modified solution for building a transition matrix from lots of chains</strong>
  <br>
  The best answer to the above question specifically was DSM's.  However, for anyone who wants to populate a transition matrix based on a list of millions of markov chains, the quickest way is this:</p>
  
  <pre><code>def fast_increment_transition_counts_from_chain(markov_chain, transition_counts_matrix):
      flat_coords = numpy.ravel_multi_index((markov_chain[:-1], markov_chain[1:]), transition_counts_matrix.shape)
      transition_counts_matrix.flat += numpy.bincount(flat_coords, minlength=transition_counts_matrix.size)
  
  def get_fake_transitions(markov_chains):
      fake_transitions = []
      for i in xrange(1,len(markov_chains)):
          old_chain = markov_chains[i - 1]
          new_chain = markov_chains[i]
          end_of_old = old_chain[-1]
          beginning_of_new = new_chain[0]
          fake_transitions.append((end_of_old, beginning_of_new))
      return fake_transitions
  
  def decrement_fake_transitions(fake_transitions, counts_matrix):
      for old_state, new_state in fake_transitions:
          counts_matrix[old_state, new_state] -= 1
  
  def fast_get_transition_counts_matrix(markov_chains, number_of_states):
      """50% faster than original, but must store 2 additional slice copies of all markov chains in memory at once.
      You might need to break up the chains into manageable chunks that don't exceed your memory.
      """
      transition_counts_matrix = numpy.zeros([number_of_states, number_of_states])
      fake_transitions = get_fake_transitions(markov_chains)
      markov_chains = list(itertools.chain(*markov_chains))
      fast_increment_transition_counts_from_chain(markov_chains, transition_counts_matrix)
      decrement_fake_transitions(fake_transitions, transition_counts_matrix)
      return transition_counts_matrix
  </code></pre>