We already discussed insertion and quick sorting in previous posts. Today, we are going to implement a quick function to perform heap sort in Python in order to answer one of Rosalind problems. Like merge sort, the worst case running time of heap sort is $\mathcal{O}(n\log n)$ The idea behind this approach is to convert an array of values into a heap, which is a (perfectly balanced) binary tree where each node is greater than each of its children and all leaves are in the leftmost position available. In the case of a heap, the children of an element n are at index 2n+1 for the left child and 2n+2 for the right child. Note that the Heap queue algorithm is readily available in Python 3, and according to the documentation, “a heapsort can be implemented by pushing all values onto a heap and then popping off the smallest values one at a time.”
Here is a recursive formulation, following the steps suggested in Cormen, Leiserson, Rivest & Stein’s Introduction to Algorithms (§6.2) in order to maintain the max-heap property (“max-heapify”):
def heapify(lst, k, n):
curr = k
left = 2 * k + 1
right = 2 * k + 2
if left < n and lst[left] > lst[curr]:
curr = left
if right < n and lst[right] > lst[curr]:
curr = right
if curr != k:
lst[curr], lst[k] = lst[k], lst[curr]
heapify(lst, curr, n)
def heapsort(lst):
n = len(lst)
for i in range(int(n/2) - 1, -1, -1):
heapify(lst, i, n)
for i in range(n - 1, 0, -1):
lst[0], lst[i] = lst[i], lst[0]
heapify(lst, 0, i)
return lst
And here is the “built-in” procedure:
from heapq import heappush, heappop
def heapsort2(iterable):
h = []
for value in iterable:
heappush(h, value)
return [heappop(h) for i in range(len(h))]
How do the two implementations compare?
In [10]: import random
In [11]: xs = [random.random() for _ in range(10**6)]
In [12]: %time heapsort(xs)
CPU times: user 14.3 s, sys: 39.7 ms, total: 14.3 s
Wall time: 14.4 s
In [13]: %time heapsort2(xs)
CPU times: user 1.4 s, sys: 21.3 ms, total: 1.42 s
Wall time: 1.43 s
So, this is clearly a big win for the later, and it suggest that we would be better using heapq for solving Rosalind problems involving heap sort. As confirmed by a bit of profiling, most of the time is spent calling heapify multiple times in our own implementation,1 while Python’s heapq allows a constant access time. For an array of size 1000, there will be 1000 “pop” and 1000 “push”, at almost a zero cost.
In [17]: import cProfile
In [18]: cProfile.run('heapsort2(xs[:1000])')
2006 function calls in 0.001 seconds
Ordered by: standard name
ncalls tottime percall cumtime percall filename:lineno(function)
1 0.000 0.000 0.001 0.001 <ipython-input-10-4d33d18db4ca>:2(heapsort2)
1 0.000 0.000 0.001 0.001 <ipython-input-10-4d33d18db4ca>:6(<listcomp>)
1 0.000 0.000 0.001 0.001 <string>:1(<module>)
1000 0.000 0.000 0.000 0.000 {built-in method _heapq.heappop}
1000 0.000 0.000 0.000 0.000 {built-in method _heapq.heappush}
1 0.000 0.000 0.001 0.001 {built-in method builtins.exec}
1 0.000 0.000 0.000 0.000 {built-in method builtins.len}
1 0.000 0.000 0.000 0.000 {method 'disable' of '_lsprof.Profiler' objects}