Implement LSD radix sort in python

Question

I'm trying to implement LSD radix sort in python, but my code used 4 loop, is there any way I can do it in three or less?

Here is the 4 loops:

First loop iterate through the input list, convert it to desired radix, and change each number to a reversed list, like so:

123 --> 173 --> [3, 7, 1] #radix == 8

Second loop iterate through the newly created list, make them have the same length by adding 0, add a pointer to the corresponding element in input list, then reverse the list.Like this:

input: [1234, 56, 7]
After first loop: [[4, 3, 2, 1], [6, 5], [7]]
After the second loop: [[0, 1, 2, 3, 4], [1, 0, 0, 5, 6], [2, 0, 0, 0, 7]]

Third loop sort element in LSD order using counting sort, code is below.

Fourth loop put every thing in a new list with sorted order using the pointers.

Here's the code for radix sort:

def radix_sort(collection, radix):

    new_collect = []
    max_len = 0

    #loop through collection, convert them to desired form and put in new_collect
    for i in collection:
        num = []
        while i >= radix:
            remain = i % radix
            num.append(remain)
            i = i // radix
        num.append(i)
        new_collect.append(num)

        if len(num) > max_len:
            max_len = len(num)

    #loop through new_collect, extend the each item to same length
    for i in range(0, len(new_collect)):
        space = max_len - len(new_collect[i])
        patch = [0 for j in range(space)]
        new_collect[i].extend(patch)

        #add a pointer to the corresponding element in collection
        new_collect[i].append(i)
        new_collect[i].reverse()

    #sort by digit with counting_sort
    for i in range(-1, -1 - max_len, -1):
        new_collect = counting_sort(new_collect, radix - 1, i)

    #create a new list with same length as collection
    return_list = list(range(len(collection)))

    #put elements in collection to return_list, using the pointer in new_collect
    for i in range(0, len(collection)):
        return_list[i] = collection[new_collect[i][0]]

    return return_list

And here is the code for counting sort:

def counting_sort(collection, d_range, digit = -1):

    d_range += 1

    B = list(range(len(collection) + 1))
    C = list(range(d_range))

    for i in range(d_range):
        C[i] = 0

    for j in collection:
        #C[j] = |{key = j}|
        C[j[digit]] += 1

    for i in range(1, d_range):
        #C[i] = |{key <= i}|
        C[i] = C[i] + C[i - 1]

    for i in range(len(collection) - 1, -1, -1):
        B[C[collection[i][digit]]] = collection[i]
        C[collection[i][digit]] -= 1

    return B[1:]

Reply to Gareth Rees 's answer.

Wow, thank you for the detailed answer!

1.I guess I'm pretty bad at naming and documenting.

2.Never even tought of using |{key < k}| that's pretty neat!

3.The second argument is not the number of keys, at least I don't have the intention to make it so, it should be the range for the input, and I add 1 so that it includes 0. So I guess it's all possible value of keys. In your improved version of the code, the value of the seq could be greater than the length of the bound, which I think will cause a out of range error.

4.About what you said in 4.1, since my radix_sort function require a radix for sorting, and I use the range of value instead of number of keys, wouldn't it be quicker to just input radix-1 or radix instead of using max()?

5.About what you said in 4.2, in my radix sort function, i need to change the raidx of input to a given value, so I need to use a loop to change it, and since I don't want the letters get involved for radixs greater than 10, I need to separate each digit, that's why I used lists. So unless I'm missing something, I don't think I can use you method in this particular scenario.

6.Do you use timeit to calculate the runtime? Or some other modules?

Thanks again for the answer!

PS: Is there anything you can say about my radix sort function's general logic? Is there anyway I can improve it?

Answer 1

1. Introduction

In this answer I'm just going to look at the function counting_sort. You'll see that there's plenty here for one answer.

2. Review

1. There's no docstring. What does this function do? What arguments does it take? What does it return?

2. When a function returns a result, it's conventional to name it with a noun or noun phrase describing the result. Hence the built-in function is called sorted rather than sort.

3. The first argument is named collection but this needs to be a sequence, so a name like seq would be better.

4. If the second argument were the number of keys (rather than the maximum key) then there would be no need to increment it (and similarly no need to decrement this argument in the callers).

5. The variables B and C could have more descriptive names, for example result and count respectively.

3. Performance

1. The lists B and C can be initialized more quickly using the * operator:

   B = [None] * (len(collection) + 1)
   C = [0] * d_range
   

2. The main loop looks like this:

   for i in range(len(collection) - 1, -1, -1):
       B[C[collection[i][digit]]] = collection[i]
       C[collection[i][digit]] -= 1
   

   But this means that collection[i] has to be looked up three times on each iteration. It is faster to iterate over the items in the collection (using reversed to do so in reverse order), avoiding the lookups:

   for item in reversed(collection):
       B[C[item[digit]]] = item
       C[item[digit]] -= 1
   

3. Now, item[digit] has to be looked up twice on each iteration. It is faster to cache this in a local variable:

   for item in reversed(collection):
       k = item[digit] # the key for this item
       B[C[k]] = item
       C[k] -= 1
   

4. Now, item[digit] is looked up twice for each item, once when counting the number of occurrences of each key, and once in the main loop. This can be reduced to a single lookup by creating an array of keys at the start:

   keys = [item[digit] for item in collection]
   

   and then using this when counting keys and in the main loop.

5. If C[k] is changed so that it has the value \$\left|\{\mathrm{key}<k\}\right|\$ rather than \$\left|\{\mathrm{key}≤k\}\right|\$, then this has three benefits: (i) the main loop can iterate over the collection in forward order, avoiding the call to reversed; (ii) the slice B[1:] at the end is avoided; (iii) the initialization of B is simplified.

Applying all these changes results in the following code:

def counting_sorted(seq, nkeys, digit=-1):
    """Return a new list containing all the items from seq in ascending
    order of item[digit], which must be a non-negative integer less
    than the parameter nkeys.

    """
    keys = [item[digit] for item in seq]

    # count[k] = |{key = k}|
    count = [0] * nkeys
    for k in keys:
        count[k] += 1

    # index[k] = |{key < k}|
    index = [0] * nkeys
    for k in range(nkeys - 1):
        index[k + 1] = index[k] + count[k]

    result = [None] * len(seq)
    for k, item in zip(keys, seq):
        result[index[k]] = item
        index[k] += 1
    return result

I find that this reduces the runtime to about 60% of the original. On my test data this is a little faster than the built-in sorted.

4. Interface improvements

There are a couple of infelicities in the interface to counting_sorted, and I think that for ease of use and understandability of the code it would be worth fixing them, even at the expense of a few percent slowdown.

1. The function requires us to pass in the value of nkeys. This is fragile (it would be easy to pass in the wrong value) and redundant (it can be computed from the data). So I would remove this parameter and add the line:

   nkeys = max(keys) + 1
   

   after the computation of keys.

2. The function looks up the sort key for each item by indexing the item with digit. But this only works for some kinds of item. Other sorting functions (like the built-in sorted) take a general key-derivation function, and I think it would be better to do so here. In your radix sort function you can make the key functions you need using operator.itemgetter

Revised code:

def counting_sorted(seq, key=None):
    """Return a new list containing all the items from seq in ascending
    order. A custom key function can be supplied to customise the sort
    order. Sort keys must be non-negative integers.

    """
    if key is None:
        keys = seq
    else:
        keys = list(map(key, seq))
    nkeys = max(keys) + 1

    # count[k] = |{key = k}|
    count = [0] * nkeys
    for k in keys:
        count[k] += 1

    # index[k] = |{key < k}|
    index = [0] * nkeys
    for k in range(nkeys - 1):
        index[k + 1] = index[k] + count[k]

    result = [None] * len(seq)
    for k, item in zip(keys, seq):
        result[index[k]] = item
        index[k] += 1
    return result

5. Answers to questions

1. Naming is one of the hardest things about programming. Reading other peoples' code is good practice—you'll find yourself continually thinking, "what on earth did the programmer mean by foo?" and then you can apply the same thought process to your own code.

2. I wanted to avoid the slice on B, and this change followed straightforwardly.

3. I think you're right and nkeys was a poor choice of name. What I meant by this was "number of (possible) distinct key values" but the name nkeys doesn't make that clear. (Naming is hard!)

4. You're right in this particular case that passing radix would be faster than calling max. But you'll see in the introduction to §4 that I gave my reasons for making these changes as "ease of use and understandability". The point is that by generalizing the function (so that it works in the general case, and not just in one specific context) it is easier to re-use, and easier for other programmers to understand and maintain.

   This is a theoretical concern—after all, maybe you'll never want to re-use this code and no other programmers will ever work on it. But considered as an exercise I think it's valuable: in a professional context it's a good idea to make code as clear and easy to maintain as possible.

5. It's easy to change radix_sort to use the revised interface. There are two sensible approaches. (i) Transform each item into a list as before, and then pass key=itemgetter(-1) on the first iteration, key=itemgetter(-2) on the second iteration, and so on. (ii) Leave the items unchanged, and pass key=lambda x: x % radix on the first iteration, key=lambda x: (x // radix) % radix on the second iteration, and so on.

6. I use timeit.

7. As I said in the introduction, this answer is quite long enough already. Maybe someone else can have a look at radix_sort?

Answer 2

I didn't find anyway to change the general logic, but with @Gareth Rees 's answer, I was able to modify the code and make the radix_sorted run about 30% faster.

1. I wasn't able to do anything to the first loop.
2. In the second loop, I avoided some repetitive look up. Removed the reverse at the end, since changing the order of running counting_sorted can do that. Also, instead of adding pointers, I added the original value into the new seq. This seems to be better in the test I did, but I'm not sure how it is gonna perform in sequence with greater numbers.
3. In the third loop, I reversed the order of sorting, so that the reverse in the second loop can be avoided.
4. The fourth loop is now merged with the return.

Here is the new code:

def radix_sorted(seq, radix):
    """
    Return a new list containing all elements from seq in ascending order.
    Second argument is required to convert seq to desired form.

    """

    converted_seq = []
    max_len = 0

    #loop through seq, convert them to desired form and put in converted_seq
    for i in seq:
        num = []
        while i >= radix:
            remain = i % radix
            num.append(remain)
            i = i // radix
        num.append(i)
        converted_seq.append(num)

        #find the max length for converted number
        if len(num) > max_len:
            max_len = len(num)

    #make every element in converted seq have same length by adding 0s
    #put original key in converted_seq
    for key, item in zip(seq, converted_seq):
        space = max_len - len(item)
        patch = [0] * space
        item.extend(patch)

        #add the element in seq to the new list
        item.append(key)

    #sort by digit with counting_sort
    for i in range(max_len):
        converted_seq = counting_sorted(converted_seq, radix - 1, i)

    #return the sorted list, using the pointer in converted_seq
    return [item[-1] for item in converted_seq]