def gap_count(pileup):
    """Given a pileup returns two matrices, one for the number of gaps and
    one for the number of non-gap reads. The first index of the matrix refers to
    forward or backward reads."""
    ngaps = pileup[:, 4, :]
    nnucl = pileup.sum(axis=1) - ngaps
    return ngaps, nnucl