def METHOD_NAME(packages_path, params_info, likelihoods_info, random_state,
                        cached=True) -> Optional[dict]:
    """
    Actual covmat finder used by `get_best_covmat`. Call directly for more control on
    the parameters used.
    Returns the same dict as `get_best_covmat`, except for the covariance matrix itself.
    """
    global _loaded_covmats_database
    covmats_database = (
            _loaded_covmats_database or get_covmat_database(packages_path, cached=cached))
    _loaded_covmats_database = covmats_database
    # Prepare params and likes aliases
    params_renames = set(chain(*[
        [p] + str_to_list(info.get("renames", [])) for p, info in
        params_info.items()]))
    likes_renames = set(chain(*[[like] + str_to_list((info or {}).get("aliases", []))
                                for like, info in likelihoods_info.items()]))
    delimiters = r"[_\.]"
    likes_regexps = [re.compile(delimiters + re.escape(_like) + delimiters)
                     for _like in likes_renames]
    # Match number of params
    score_params = (
        lambda covmat: len(set(covmat["params"]).intersection(params_renames)))
    best_p = get_best_score(covmats_database, score_params)
    if not best_p:
        log.warning(
            "No covariance matrix found including at least one of the given parameters")
        return None
    # Match likelihood names / keywords
    # No debug print here: way too many!
    score_likes = (
        lambda covmat: len([0 for r in likes_regexps if r.search(covmat["name"])]))
    best_p_l = get_best_score(best_p, score_likes)
    if is_debug(log):
        log.debug("Subset based on params + likes:\n - " +
                  "\n - ".join([b["name"] for b in best_p_l]))
    # Finally, in case there is more than one, select shortest #params and name (simpler!)
    # #params first, to avoid extended models with shorter covmat name
    def score_simpler_params(covmat):
        return -len(covmat["params"])
    best_p_l_sp = get_best_score(best_p_l, score_simpler_params)
    if is_debug(log):
        log.debug("Subset based on params + likes + fewest params:\n - " +
                  "\n - ".join([b["name"] for b in best_p_l_sp]))
    score_simpler_name = (
        lambda covmat: -len(covmat["name"].replace("_", " ").replace("-", " ").split()))
    best_p_l_sp_sn = get_best_score(best_p_l_sp, score_simpler_name)
    if is_debug(log):
        log.debug("Subset based on params + likes + fewest params + shortest name:\n - " +
                  "\n - ".join([b["name"] for b in best_p_l_sp_sn]))
    # if there is more than one (unlikely), just pick one at random
    if len(best_p_l_sp_sn) > 1:
        log.warning("WARNING: >1 possible best covmats: %r",
                    [b["name"] for b in best_p_l_sp_sn])
    random_state = np.random.default_rng(random_state)
    return best_p_l_sp_sn[random_state.choice(range(len(best_p_l_sp_sn)))].copy()