docs: def posthoc_quade(a, y_col=None, block_col=None, group_col=None, dist='t', melted=False, sort=False, p_adjust=None):
  
      '''Calculate pairwise comparisons using Quade's post hoc test for
      unreplicated blocked data. This test is usually conducted if significant
      results were obtained by the omnibus test [1]_, [2]_, [3]_.
  
      Parameters
      ----------
      a : array_like or pandas DataFrame object
          An array, any object exposing the array interface or a pandas
          DataFrame.
  
          If `melted` is set to False (default), `a` is a typical matrix of
          block design, i.e. rows are blocks, and columns are groups. In this
          case you do not need to specify col arguments.
  
          If `a` is an array and `melted` is set to True,
          y_col, block_col and group_col must specify the indices of columns
          containing elements of correspondary type.
  
          If `a` is a Pandas DataFrame and `melted` is set to True,
          y_col, block_col and group_col must specify columns names (string).
  
      y_col : str or int, optional
          Must be specified if `a` is a pandas DataFrame object.
          Name of the column that contains y data.
  
      block_col : str or int
          Must be specified if `a` is a pandas DataFrame object.
          Name of the column that contains blocking factor values.
  
      group_col : str or int
          Must be specified if `a` is a pandas DataFrame object.
          Name of the column that contains treatment (group) factor values.
  
      dist : str, optional
          Method for determining p values.
          The default distribution is "t", else "normal".
  
      melted : bool, optional
          Specifies if data are given as melted columns "y", "blocks", and
          "groups".
  
      sort : bool, optional
          If True, sort data by block and group columns.
  
      p_adjust : str, optional
          Method for adjusting p values.
          See statsmodels.sandbox.stats.multicomp for details.
          Available methods are:
          'bonferroni' : one-step correction
          'sidak' : one-step correction
          'holm-sidak' : step-down method using Sidak adjustments
          'holm' : step-down method using Bonferroni adjustments
          'simes-hochberg' : step-up method  (independent)
          'hommel' : closed method based on Simes tests (non-negative)
          'fdr_bh' : Benjamini/Hochberg  (non-negative)
          'fdr_by' : Benjamini/Yekutieli (negative)
          'fdr_tsbh' : two stage fdr correction (non-negative)
          'fdr_tsbky' : two stage fdr correction (non-negative)
  
      Returns
      -------
      Pandas DataFrame containing p values.
  
      References
      ----------
      .. [1] W. J. Conover (1999), Practical nonparametric Statistics, 3rd. Edition,
          Wiley.
  
      .. [2] N. A. Heckert and J. J. Filliben (2003). NIST Handbook 148: Dataplot
          Reference Manual, Volume 2: Let Subcommands and Library Functions.
          National Institute of Standards and Technology Handbook Series, June 2003.
  
      .. [3] D. Quade (1979), Using weighted rankings in the analysis of complete
          blocks with additive block effects. Journal of the American Statistical
          Association, 74, 680-683.
  
      Examples
      --------
      >>> x = np.array([[31,27,24],[31,28,31],[45,29,46],[21,18,48],[42,36,46],[32,17,40]])
      >>> sp.posthoc_quade(x)
      '''
  
      if melted and not all([block_col, group_col, y_col]):
          raise ValueError('block_col, group_col, y_col should be explicitly specified if using melted data')
  
      def compare_stats_t(i, j):
          dif = np.abs(S[groups[i]] - S[groups[j]])
          tval = dif / denom
          pval = 2. * ss.t.sf(np.abs(tval), df = (b - 1) * (k - 1))
          return pval
  
      def compare_stats_norm(i, j):
          dif = np.abs(W[groups[i]] * ff - W[groups[j]] * ff)
          zval = dif / denom
          pval = 2. * ss.norm.sf(np.abs(zval))
          return pval
  
      x, y_col, group_col, block_col = __convert_to_block_df(a, y_col, group_col, block_col, melted)
  
      if not sort:
          x[group_col] = Categorical(x[group_col], categories=x[group_col].unique(), ordered=True)
          x[block_col] = Categorical(x[block_col], categories=x[block_col].unique(), ordered=True)
      x.sort_values(by=[block_col, group_col], ascending=True, inplace=True)
      x.dropna(inplace=True)
  
      groups = x[group_col].unique()
      k = len(groups)
      b = x[block_col].unique().size
  
      x['r'] = x.groupby(block_col)[y_col].rank()
      q = (x.groupby(block_col)[y_col].max() - x.groupby(block_col)[y_col].min()).rank()
      x['rr'] = x['r'] - (k + 1)/2
      x['s'] = x.apply(lambda x, y: x['rr'] * y[x['blocks']], axis=1, args=(q,))
      x['w'] = x.apply(lambda x, y: x['r'] * y[x['blocks']], axis=1, args=(q,))
      A = (x['s'] ** 2).sum()
      S = x.groupby(group_col)['s'].sum()
      B = np.sum(S ** 2) / b
      W = x.groupby(group_col)['w'].sum()
  
      vs = np.zeros((k, k), dtype=np.float)
      combs = it.combinations(range(k), 2)
  
      tri_upper = np.triu_indices(vs.shape[0], 1)
      tri_lower = np.tril_indices(vs.shape[0], -1)
      vs[:,:] = 0
  
      if dist == 't':
          denom = np.sqrt((2 * b * (A - B)) / ((b - 1) * (k - 1)))
  
          for i, j in combs:
              vs[i, j] = compare_stats_t(i, j)
  
      else:
          n = b * k
          denom = np.sqrt((k * (k + 1) * (2 * n + 1) * (k-1)) / (18 * n * (n + 1)))
          ff = 1. / (b * (b + 1)/2)
  
          for i, j in combs:
              vs[i, j] = compare_stats_norm(i, j)
  
      if p_adjust:
          vs[tri_upper] = multipletests(vs[tri_upper], method = p_adjust)[1]
  
      vs[tri_lower] = vs.T[tri_lower]
      np.fill_diagonal(vs, -1)
      return DataFrame(vs, index=groups, columns=groups)