Validating Input, Initializing the Tree, and What Slices Copy

Code along video

Learning objectives

In this code along, we will write the subroutines of upgma() that check our input and initialize the tree. Then we will find out what really happens when we take a slice of a list, which will tell us whether our clusters list behaves the way we need it to.

Code along summary

Setup

We will continue working in the upgma directory from the previous code along, in functions.py and main.py.

Validating input

To assert that a matrix is square, we first check that it is not empty. Then we range over all of its rows and make sure that each row has length equal to the number of rows; a list of lists might have rows of different lengths, so this check is necessary. If we survive these checks, we do not need to take any action.

def assert_square_matrix(mtx: DistanceMatrix) -> None:
    """
    Validate that a distance matrix is square.

    Args:
        mtx (DistanceMatrix): The matrix to validate.

    Raises:
        ValueError: If the matrix is empty or not square.
    """
    num_rows = len(mtx)
    if num_rows == 0:
        raise ValueError("Empty matrix given.")

    # idea: range over all rows, and each row should have length == num_rows
    for i in range(num_rows):
        if len(mtx[i]) != num_rows:
            raise ValueError("Error: non-square matrix given.")

    # if we make it here, we know the matrix is square.

Checking that the number of rows of the matrix matches the number of species is even easier.

def assert_same_number_species(mtx: DistanceMatrix, species_names: list[str]) -> None:
    """
    Validate that the number of species matches the matrix dimension.

    Args:
        mtx (DistanceMatrix): Square distance matrix.
        species_names (list[str]): Species labels.

    Raises:
        ValueError: If their sizes do not match.
    """
    if len(mtx) != len(species_names):
        raise ValueError("Error: number of rows in matrix not equal to # species.")

Initializing the tree

Our initialize_tree() function takes the list of species names and returns a tree with 2n – 1 nodes, where n is the number of species. We start with an empty list, which we annotate as a Tree to remind ourselves what it represents. Then we create each node, which gets a default age of 0, and append it to the tree. If the node is one of the first num_leaves nodes, then it is a leaf, and its label is the corresponding species name; otherwise, it is an internal node, and we label it with an f-string containing its index.

def initialize_tree(species_names: list[str]) -> Tree:
    """
    Allocate all 2n - 1 nodes for a tree with n leaves.

    Args:
        species_names (list[str]): Names for the n leaf nodes.

    Returns:
        Tree: A Tree of labeled leaves followed by labeled internal nodes.
    """
    num_leaves = len(species_names)

    t: Tree = []

    # let's give the tree 2*n-1 total nodes
    for i in range(2*num_leaves - 1):
        new_node = Node()  # default age = 0
        t.append(new_node)

        # we give it a label too
        if i < num_leaves:
            t[i].label = species_names[i]
        else:
            # internal node, so we're just going to label it with i
            t[i].label = f"Ancestor Species: {i}"

    return t

We have not set any edges, and every node has age 0. The ages of internal nodes and the edges connecting nodes will be set as UPGMA runs.

What does a slice copy?

Let’s return to the line clusters = t[:num_leaves] in upgma(). Another way of asking our question from the previous lesson is: should this line make new nodes? To find out, let’s run an experiment in main(), starting with a list of integers. We set b equal to the first two elements of a, change the first element of b, and print a.

def main() -> None:
    print("Happy trees.")

    a = [10, 20, 30, 40]
    b = a[:2]
    print(b)

    b[0] = 99
    print("a is", a)
    print("b is", b)

Navigate into python/src/upgma and run python3 main.py. We see that a did not change, and so b holds a duplicated form of the data.

Happy trees.
[10, 20]
a is [10, 20, 30, 40]
b is [99, 20]

Now let’s try the same thing with a list of nodes. Since Node is declared in datatypes.py, we first need to add from datatypes import Node to the top of main.py. Then add the following lines to the end of main().

a = [Node(label="Axolotl"), Node(label="Blobfish")]
b = a[:1]
print(b[0].label)
b[0].label = "Chicken"
print("the first label of a is", a[0].label)

When we run our code, we see that the label of the first node of a has changed!

Axolotl
the first label of a is Chicken

The magic word is mutability. When we take a slice of a list of integers, Python duplicates the integers. When we take a slice of a list of nodes, Python also makes copies, but what it copies are references to the nodes. Think of the memory deep in the bowels of your computer as a storage facility. The node with the label Axolotl lives in a storage locker, and the list a holds a key to that locker. When we take the slice a[:1], Python makes a copy of the key, and with a copy of the key, b can unlock the locker and swap the axolotl for a chicken. So the behavior of slicing is the same in both cases (Python copies what is in the list), but the end result differs because a list of nodes contains references to mutable objects.

# when we take a sublist (slice) of a collection of mutable objects, we get access to those objects (can change the original list)

Now we can answer our question about upgma(). The line clusters = t[:num_leaves] does not make new nodes. It makes copies of the keys to the lockers that contain the leaves of the tree, which is exactly what we want: clusters[0] refers to the very same node as t[0], and so on. Let’s update our comment in upgma() to remind ourselves.

clusters = t[:num_leaves]  # is this OK? Should this make *new* nodes? No! It makes copies of the keys to the locker that contains the leaves

The same idea applies to the line t[p].child1 = clusters[row]. This assignment does not create a copy of a node, which would be a disaster for our tree. Instead, it hands t[p].child1 a copy of the key to the node in clusters[row].

t[p].child1 = clusters[row]  # this means update child 1 of current node to be the reference to the node in clusters[row] ... NOT "make a new node"

In many languages, this behavior is handled with something called pointers, which take a little effort to learn at first but are easier to reason about in the long run. In Python, you just need to remember that if you take a slice of a list whose elements are mutable, you get references to the original objects instead of duplicates of them.

Looking ahead

We have validated our input, initialized our tree, and confirmed that clusters refers to the leaves of the tree. In the next code along, we will write the subroutines that find the minimum element of the matrix and that shrink the matrix and the clusters at each step.