Repository: scipy/scipy
File: scipy/cluster/hierarchy.py (Python)

Code before:
    """
    Return the maximum distance between any non-singleton cluster.

    Parameters
    ----------
    Z : ndarray
        The hierarchical clustering encoded as a matrix. See
        ``linkage`` for more information.

    Returns
    -------
    maxdists : ndarray
        A ``(n-1)`` sized numpy array of doubles; ``MD[i]`` represents
        the maximum distance between any cluster (including
        singletons) below and including the node with index i. More
        specifically, ``MD[i] = Z[Q(i)-n, 2].max()`` where ``Q(i)`` is the
        set of all node indices below and including node i.

    See Also
    --------
    linkage : for a description of what a linkage matrix is.
    is_monotonic : for testing for monotonicity of a linkage matrix.

    Examples
    --------
    >>> from scipy.cluster.hierarchy import median, maxdists
    >>> from scipy.spatial.distance import pdist

    Given a linkage matrix ``Z``, `scipy.cluster.hierarchy.maxdists`
    computes for each new cluster generated (i.e., for each row of the linkage
    matrix) what is the maximum distance between any two child clusters.

    Due to the nature of hierarchical clustering, in many cases this is going
    to be just the distance between the two child clusters that were merged
    to form the current one - that is, Z[:,2].

    However, for non-monotonic cluster assignments such as
    `scipy.cluster.hierarchy.median` clustering this is not always the
    case: There may be cluster formations were the distance between the two
    clusters merged is smaller than the distance between their children.

    We can see this in an example:

    >>> X = [[0, 0], [0, 1], [1, 0],
    ...      [0, 4], [0, 3], [1, 4],
    ...      [4, 0], [3, 0], [4, 1],
    ...      [4, 4], [3, 4], [4, 3]]

    >>> Z = median(pdist(X))
    >>> Z
    array([[ 0.        ,  1.        ,  1.        ,  2.        ],

Review comment (suggestion):
tested up to 1m observations, 5m nodes. It just never seems to catch up.
As if `unique_all` was not implemented on top of a hashmap...

Also, `unique_all` carries compatibility headaches. Dask and JAX support require _vigorous massaging_, PyTorch doesn't implement it at all, and in general it is not mandatory in the Array API.

Code after:
    cluster assignments are isomorphic:

    >>> is_isomorphic(T, T_)
    True

    """
    xp = array_namespace(T1, T2)
    T1 = _asarray(T1, xp=xp)
    T2 = _asarray(T2, xp=xp)

    if T1.ndim != 1:
        raise ValueError('T1 must be one-dimensional.')
    if T2.ndim != 1:
        raise ValueError('T2 must be one-dimensional.')
    if T1.shape != T2.shape:
        raise ValueError('T1 and T2 must have the same number of elements.')

    # For each pair of (i, j) indices, test that
    # T1[i] == T1[j] <--> T2[i] == T2[j]
    # In other words, if the same symbol appears multiple times in T1,
    # then a potentially different symbol must also be repeated in the same
    # positions in T2, and vice versa.

    # O(n*log(n)) algorithm.
    # It is also possible to write a O(n) algorithm on top of unique_all(),
    # but in practice it's much slower even for large clusters.
    idx = xp.argsort(T1)
    T1 = xp.take(T1, idx)
    T2 = xp.take(T2, idx)
    changes1 = T1 != xp.roll(T1, -1)
    changes2 = T2 != xp.roll(T2, -1)
    return xp.all(changes1 == changes2)


@lazy_cython
def maxdists(Z):
    """
    Return the maximum distance between any non-singleton cluster.

    Parameters
    ----------
    Z : ndarray
        The hierarchical clustering encoded as a matrix. See
        ``linkage`` for more information.

    Returns
    -------
    maxdists : ndarray
        A ``(n-1)`` sized numpy array of doubles; ``MD[i]`` represents
        the maximum distance between any cluster (including
        singletons) below and including the node with index i. More