def posthoc_mackwolfe(a, val_col, group_col, p=None, n_perm=100, sort=False, p_adjust=None):

    '''Mack-Wolfe Test for Umbrella Alternatives.

    In dose-finding studies one may assume an increasing treatment effect with
    increasing dose level. However, the test subject may actually succumb to
    toxic effects at high doses, which leads to decresing treatment
    effects [1]_, [2]_.

    The scope of the Mack-Wolfe Test is to test for umbrella alternatives for
    either a known or unknown point P (i.e. dose-level), where the peak
    (umbrella point) is present.

    Parameters
    ----------
    a : array_like or pandas DataFrame object
        An array, any object exposing the array interface or a pandas
        DataFrame.

    val_col : str, optional
        Name of a DataFrame column that contains dependent variable values (test
        or response variable). Values should have a non-nominal scale. Must be
        specified if `a` is a pandas DataFrame object.

    group_col : str, optional
        Name of a DataFrame column that contains independent variable values
        (grouping or predictor variable). Values should have a nominal scale
        (categorical). Must be specified if `a` is a pandas DataFrame object.

    p : int, optional
        The a-priori known peak as an ordinal number of the treatment group
        including the zero dose level, i.e. p = {1, ..., k}. Defaults to None.

    sort : bool, optional
        If True, sort data by block and group columns.

    Returns
    -------
    p : float
        P value.
    stat : float
        Statistic.

    References
    ----------
    .. [1] Chen, I.Y. (1991) Notes on the Mack-Wolfe and Chen-Wolfe Tests for
        Umbrella Alternatives. Biom. J., 33, 281-290.
    .. [2] Mack, G.A., Wolfe, D. A. (1981) K-sample rank tests for umbrella
        alternatives. J. Amer. Statist. Assoc., 76, 175-181.

    Examples
    --------
    >>> x = np.array([[10,'a'], [59,'a'], [76,'b'], [10, 'b']])
    >>> sp.posthoc_mackwolfe(x, val_col = 0, group_col = 1)

    '''

    x = __convert_to_df(a, val_col, group_col)

    if not sort:
        x[group_col] = Categorical(x[group_col], categories=x[group_col].unique(), ordered=True)
    x.sort_values(by=[group_col], ascending=True, inplace=True)

    k = x[group_col].unique().size

    if p:
        if p > k:
            print("Selected 'p' > number of groups:", str(p), " > ", str(k))
            return False
        elif p < 1:
            print("Selected 'p' < 1: ", str(p))
            return False

    Rij = x[val_col].rank()
    n = x.groupby(group_col)[val_col].count()

    def _fn(Ri, Rj):
        return np.sum(Ri.apply(lambda x: Rj[Rj > x].size))

    def _ustat(Rij, g, k):
        levels = np.unique(g)
        U = np.identity(k)

        for i in range(k):
            for j in range(i):
                U[i,j] = _fn(Rij[x[group_col] == levels[i]], Rij[x[group_col] == levels[j]])
                U[j,i] = _fn(Rij[x[group_col] == levels[j]], Rij[x[group_col] == levels[i]])

        return U

    def _ap(p, U):
        tmp1 = 0
        if p > 0:
            for i in range(p):
                for j in range(i+1, p+1):
                    tmp1 += U[i,j]
        tmp2 = 0
        if p < k:
            for i in range(p, k):
                for j in range(i+1, k):
                    tmp2 += U[j,i]

        return tmp1 + tmp2

    def _n1(p, n):
        return np.sum(n[:p+1])

    def _n2(p, n):
        return np.sum(n[p:k])

    def _mean_at(p, n):
        N1 = _n1(p, n)
        N2 = _n2(p, n)
        return (N1**2 + N2**2 - np.sum(n**2) - n.iloc[p]**2)/4

    def _var_at(p, n):
        N1 = _n1(p, n)
        N2 = _n2(p, n)
        N = np.sum(n)

        var = (2 * (N1**3 + N2**3) + 3 * (N1**2 + N2**2) -\
                np.sum(n**2 * (2*n + 3)) - n.iloc[p]**2 * (2 * n.iloc[p] + 3) +\
                12. * n.iloc[p] * N1 * N2 - 12. * n.iloc[p] ** 2 * N) / 72.
        return var

    if p:
        if (x.groupby(val_col).count() > 1).any().any():
            print("Ties are present")
        U = _ustat(Rij, x[group_col], k)
        est = _ap(p, U)
        mean = _mean_at(p, n)
        sd = np.sqrt(_var_at(p, n))
        stat = (est - mean)/sd
        p_value = ss.norm.sf(stat)
    else:
        U = _ustat(Rij, x[group_col], k)
        Ap = np.array([_ap(i, U) for i in range(k)]).ravel()
        mean = np.array([_mean_at(i, n) for i in range(k)]).ravel()
        var = np.array([_var_at(i, n) for i in range(k)]).ravel()
        A = (Ap - mean) / np.sqrt(var)
        stat = np.max(A)
        p = A == stat
        est = None

        mt = []
        for i in range(n_perm):

            ix = Series(np.random.permutation(Rij))
            Uix = _ustat(ix, x[group_col], k)
            Apix = np.array([_ap(i, Uix) for i in range(k)])
            Astarix = (Apix - mean) / np.sqrt(var)
            mt.append(np.max(Astarix))

        mt = np.array(mt)
        p_value = mt[mt > stat] / n_perm

    return p_value, stat