Returns an array whose ith row is the q[i]th quantile of the
    function.

    :Arguments:
        func_stacks: The samples of the function. func_stacks[i,:]
            gives sample i.
        qlist: A list or array of the quantiles you would like.

    :SeeAlso: func_envelopes, func_hist, weightplot

def func_quantiles(node, qlist=(.025, .25, .5, .75, .975)):
    """
    Returns an array whose ith row is the q[i]th quantile of the
    function.

    :Arguments:
        func_stacks: The samples of the function. func_stacks[i,:]
            gives sample i.
        qlist: A list or array of the quantiles you would like.

    :SeeAlso: func_envelopes, func_hist, weightplot
    """

    # For very large objects, this will be rather long.
    # Too get the length of the table, use obj.trace.length()

    if isinstance(node, Variable):
        func_stacks = node.trace()
    else:
        func_stacks = node

    if any(qlist < 0.) or any(qlist > 1.):
        raise TypeError('The elements of qlist must be between 0 and 1')

    func_stacks = func_stacks.copy()

    N_samp = shape(func_stacks)[0]
    func_len = tuple(shape(func_stacks)[1:])

    func_stacks.sort(axis=0)

    quants = zeros((len(qlist), func_len), dtype=float)
    alphas = 1. - abs(array(qlist) - .5) / .5

    for i in range(len(qlist)):
        quants[i, ] = func_stacks[int(qlist[i] * N_samp), ]

    return quants, alphas