def probvec(m, k, random_state=None, parallel=True):
    """
    Return m randomly sampled probability vectors of dimension k.

    Parameters
    ----------
    m : scalar(int)
        Number of probability vectors.

    k : scalar(int)
        Dimension of each probability vectors.

    random_state : int or np.random.RandomState, optional
        Random seed (integer) or np.random.RandomState instance to set
        the initial state of the random number generator for
        reproducibility. If None, a randomly initialized RandomState is
        used.

    parallel : bool(default=True)
        Whether to use multi-core CPU (parallel=True) or single-threaded
        CPU (parallel=False). (Internally the code is executed through
        Numba.guvectorize.)

    Returns
    -------
    x : ndarray(float, ndim=2)
        Array of shape (m, k) containing probability vectors as rows.

    Examples
    --------
    >>> qe.random.probvec(2, 3, random_state=1234)
    array([[ 0.19151945,  0.43058932,  0.37789123],
           [ 0.43772774,  0.34763084,  0.21464142]])

    """
    if k == 1:
        return np.ones((m, k))

    # if k >= 2
    random_state = check_random_state(random_state)
    r = random_state.random_sample(size=(m, k-1))
    x = np.empty((m, k))

    # Parse Parallel Option #
    if parallel:
        _probvec_parallel(r, x)
    else:
        _probvec_cpu(r, x)

    return x