'Compute quantiles for a weighted sample.
Parameters
probs : array-like
A vector of probability points at which to calculate the
quantiles.  Each element of `probs` should fall in [0, 1].
return_pandas : bool
If True, return value is a Pandas DataFrame or Series.
Otherwise returns a ndarray.
Returns
quantiles : Series, DataFrame, or ndarray
If `return_pandas` = True, returns one of the following:
* data are 1d, `return_pandas` = True: a Series indexed by
the probability points.
* data are 2d, `return_pandas` = True: a DataFrame with
the probability points as row index and the variables
as column index.
If `return_pandas` = False, returns an ndarray containing the
same values as the Series/DataFrame.
Notes
To compute the quantiles, first, the weights are summed over
exact ties yielding distinct data values y_1 < y_2 < ..., and
corresponding weights w_1, w_2, ....  Let s_j denote the sum
of the first j weights, and let W denote the sum of all the
weights.  For a probability point p, if pW falls strictly
between s_j and s_{j+1} then the estimated quantile is
y_{j+1}.  If pW = s_j then the estimated quantile is (y_j +
y_{j+1})/2.  If pW < p_1 then the estimated quantile is y_1.
References
SAS documentation for weighted quantiles:
https://support.sas.com/documentation/cdl/en/procstat/63104/HTML/default/viewer.htm#procstat_univariate_sect028.htm'

def quantile(self, probs, return_pandas=True):

    import pandas as pd
    probs = np.asarray(probs)
    probs = np.atleast_1d(probs)
    if (self.data.ndim == 1):
        rslt = self._quantile(self.data, probs)
        if return_pandas:
            rslt = pd.Series(rslt, index=probs)
    else:
        rslt = []
        for vec in self.data.T:
            rslt.append(self._quantile(vec, probs))
        rslt = np.column_stack(rslt)
        if return_pandas:
            columns = [('col%d' % (j + 1)) for j in range(rslt.shape[1])]
            rslt = pd.DataFrame(data=rslt, columns=columns, index=probs)
    if return_pandas:
        rslt.index.name = 'p'
    return rslt
