def max_variance_genes(data, nbins=5, frac=0.2):
    indices = []
    if sparse.issparse(data):
        means, var = sparse_mean_var(data)
    else:
        means = data.mean(1)
        var = data.var(1)
    mean_indices = means.argsort()
    n_elements = int(data.shape[0]/nbins)
    frac_elements = int(n_elements*frac)
    for i in range(nbins):
        bin_i = mean_indices[i*n_elements : (i+1)*n_elements]
        if i==nbins-1:
            bin_i = mean_indices[i*n_elements :]
        var_i = var[bin_i]
        var_sorted = var_i.argsort()
        top_var_indices = var_sorted[len(bin_i) - frac_elements:]
        ind = bin_i[top_var_indices]
        ind = [index for index in ind if var[index]>0]
        indices.extend(ind)
    return indices