Performances et optimisation de métaheuristique
Workshop

Guide tuteur

Version	Date	Concepteur	Relecteur	Commentaires
1.1	2022-09-27	Benjamin COHEN BOULAKIA		Réorganisation du Workshop sous forme de points de discussion avec les étudiants Quelques explications plus détaillées
1.0	2022-03-15	Benjamin COHEN BOULAKIA	Jean-François DOLLINGER	Version initiale, reprenant la partie Qualité de la solution du Workshop Métaheuristiques et la complétant avec l'étude du comportement statistique de l'algorithme (tabou simple).

Compétences du Workshop

• Expérimenter une ou plusieurs bibliothèques de calcul scientifique
• Valider et Comparer les résultats expérimentaux obtenus par différentes approches
• Recommander des optimisations des algorithmes implémentés

Objectifs

Ce Workshop reprend le problème du sac à dos unidimensionnel en 0–1 étudié dans le Workshop précédent, ainsi que la méthode par recherche avec liste tabou implémentée pour y répondre. Les étudiants vont examiner expérimentalement le comportement de la recherche tabou, interpréter les résultats obtenus, et identifier des opportunités d'optimisation sur cette base. Les étudiants doivent éditer ce Notebook en répondant aux questions posées, et en complétant le code fourni. Ils peuvent aussi exporter ensuite le code Python pour exécution hors de Jupyter.

Par rapport au Workshop précédent, le code à fournir est nettement moins compliqué. L'algorithme est déjà implémenté, il s'agira essentiellement de coder des tests de fonctionnement pour examiner le comportement de cet algorithme. Aucune difficulté particulière en terme de développement.

Philosophie et recommandations

Le but de ce Workshop est d'amener les étudiants à adopter un regard critique concernant leur algorithme, à déterminer la démarche analytique pour interpréter les résultats expérimentaux, et de prendre du recul par rapport à ces méthodes. Les dernières étapes les amènent à travailler sur des scénarios de test simples, et à réfléchir à des scénarios plus élaborés. Il faut les encourager à mener ce travail complémentaire en autonomie, car c'est ce type d'étude statistique qui sera attendu d'eux pour leur projet.

Certaines expériences prennent du temps. Ce Workshop s'appuie donc sur la version optimisée de la recherche tabou du Workshop précédent, qu'il est nécessaire d'avoir implémentée pour obtenir des résultats pertinents. Des CPU de type Core (Intel) ou Ryzen (AMD) sont recommandés, les solutions basse consommation (type Celeron) risquant de poser des problèmes de performance. Le cas échéant, il faudra réduire les paramètres de l'expérience (notamment le nombre de tests), quitte à obtenir des résultats non représentatifs. Il pourra aussi être utile d'afficher sur un écran les courbes de la correction pour en discuter.

Par ailleurs, de nombreux liens sont égrenés tout au long du Workshop, il faut encourager les étudiants à les visiter systématiquement.

1. Introduction

Nous reprenons pour la dernière fois le problème du sac à dos, et on va chercher à savoir si les algorithmes qu'on a conçus sont efficaces

Explications On a implémenté la dernière fois trois métaheuristiques pour traiter ce problème, qui construisent des solutions valides à une instance du problème. Mais on ne sait pas ce que valent ces solutions qu'on construit. Sont-elles de bonne qualité ? Peut-on faire mieux ?

Et puis, on ne sais pas non plus si ce qu'on obtient est représentatif. Après tout, on n'a traité qu'une seule instance. Peut-être qu'on a eu de la chance, et que cette instance était très facile. Ou au contraire, peut-être que les solutions qu'on a construit étaient mauvaises, mais qu'on n'a pas eu de chance et qu'en général nos algorithmes se comportent plutôt bien.

On va donc mettre en place un plan d'expérience pour répondre à ces questions. Tout ce workshop consiste à réfléchir à ce que doit faire ce plan d'expérience.

Ce qu'il faudrait, c'est générer des instances aléatoires pour pouvoir produire quelques statistiques descriptives qui nous donneraient une idée de tout ça. On pourrait faire varier certains paramètres des instances générées (nombre d'objets, taille ou valeur maximale…) pour voir leur impact sur les algorithmes

Et on peut aussi faire varier les paramètres des algorithmes eux-mêmes. Sur la recherche simple (sans multi-start) par exemple, c'est la taille de la liste tabou et le nombre d'itérations qui sont les principaux moyens d'améliorer la qualité des solutions calculées. Justement, on va rester sur cette recherche tabou simple, sans faire de multi-start. À votre avis, pourquoi ?

Le multi-start est lui aussi intéressant à tester, mais si on mélange les deux, on aura du mal à interpréter les résultats. Il faut d'abord étudier plus en profondeur le comportement de la recherche tabou en elle-même. Ensuite, une fois qu'on aura une vision plus claire de son comportement, on pourra l'intégrer dans une étude du mutl-start.

En plus, faire du multi-start rallonge considérablement les calculs.

Dans ce plan d'expérience, il va être très important que vos résultats soient reproductibles, si on veut pouvoir comparer les résultats qu'on obtient. Imaginez qu'on pense avoir amélioré notre algorithme en modifiant de ses paramètres, alors qu'il est juste tombé sur une instance plus facile à résoudre.

C'est la raison pour laquelle on utilise pour notre génération aléatoire une graine explicite. On réutilisera cette valeur pour chaque expérience avec génération aléatoire, que ce soit une génération d'instance ou de solution initiale (mais pas pour chaque essai dans notre expérience, évidemment).

Commencez par intégrer le code de la recherche tabou que vous avez conçue la dernière fois.

from collections import deque
import random
from functools import lru_cache

def random_objets(poids_max, val_max):
    """
    Cette fonction génère des objets de poids et de valeur
    aléatoires (bornés par les valeurs passées en parametre).
    
    Renvoie un tuple de 2 dictionnaires (poids,valeur)
    """
    poids_objets =\
        {i: random.randint(1, poids_max) for i in range(nb_objets)} #SOLUTION
    valeur_objets =\
        {i: random.randint(1, val_max) for i in range(nb_objets)}   #SOLUTION
                                                                    #SOLUTION
    return poids_objets, valeur_objets                              #SOLUTION

def random_solution():
    sac = tuple(random.choice([False, True]) for _ in range(nb_objets)) #SOLUTION
                                                                        #SOLUTION
    while (poids_contenu(sac) > capacite):                              #SOLUTION
        objets_presents = tuple(                                        #SOLUTION
            i for i, val in enumerate(sac) if val)                      #SOLUTION
        objet_supprime = random.choice(objets_presents)                 #SOLUTION
        sac = sac[:objet_supprime] + (False,) + sac[objet_supprime+1:]  #SOLUTION
                                                                        #SOLUTION
    return sac                                                          #SOLUTION

@lru_cache(maxsize=None)
def poids_contenu(sac):
    """
    Cette fonction renvoie la somme des poids des objets dans le sac
    """
    return sum(poids_objets[i] for i, val in enumerate(sac) if val) #SOLUTION

@lru_cache(maxsize=None)
def valeur_contenu(sac):
    """
    Cette fonction renvoie la somme des valeurs des objets dans le sac
    """
    return sum(valeur_objets[i] for i, val in enumerate(sac) if val) #SOLUTION

@lru_cache(maxsize=None)
def voisinage(sac):
    voisins = []                                                     #SOLUTION
    for k in range(len(sac)):                                        #SOLUTION
        voisin = sac[:k] + (not(sac[k]),) + sac[k+1:]                #SOLUTION
        if (poids_contenu(voisin) <= capacite):                      #SOLUTION
            voisins.append(voisin)                                   #SOLUTION
    return voisins                                                   #SOLUTION

def recherche_tabou(solution_initiale, taille_tabou, iter_max):
    """
    1. On part d'un élément de notre ensemble de recherche qu'on déclare élément courant
    2. On considère le voisinage de l'element courant et on choisit le  meilleur d'entre
       eux comme nouvel element courant, parmi ceux absents de la liste tabou, et on l'ajoute
       a la liste tabou
    3. On boucle jusqu'a condition de sortie.
    """
    nb_iter = 0                                                                #SOLUTION
    liste_tabou = deque((), maxlen = taille_tabou)                             #SOLUTION
                                                                               #SOLUTION
    # variables solutions pour la recherche du voisin optimal non tabou        #SOLUTION
    solution_courante = solution_initiale                                      #SOLUTION
    meilleure = solution_initiale                                              #SOLUTION
    meilleure_globale = solution_initiale                                      #SOLUTION
                                                                               #SOLUTION
    # variables valeurs pour la recherche du voisin optimal non tabou          #SOLUTION
    valeur_meilleure = valeur_contenu(solution_initiale)                       #SOLUTION
    valeur_meilleure_globale = valeur_meilleure                                #SOLUTION
                                                                               #SOLUTION
    while (nb_iter < iter_max):                                                #SOLUTION
        valeur_meilleure = -1                                                  #SOLUTION
                                                                               #SOLUTION
        # on parcourt tous les voisins de la solution courante                 #SOLUTION
        for voisin in voisinage(solution_courante):                            #SOLUTION
            valeur_voisin=valeur_contenu(voisin)                               #SOLUTION
                                                                               #SOLUTION
            # MaJ meilleure solution non taboue trouvée                        #SOLUTION
            if valeur_voisin > valeur_meilleure and voisin not in liste_tabou: #SOLUTION
                valeur_meilleure = valeur_voisin                               #SOLUTION
                meilleure = voisin                                             #SOLUTION
                                                                               #SOLUTION
        # on met à jour la meilleure solution rencontrée depuis le début       #SOLUTION
        if valeur_meilleure > valeur_meilleure_globale:                        #SOLUTION
            meilleure_globale = meilleure                                      #SOLUTION
            valeur_meilleure_globale = valeur_meilleure                        #SOLUTION
            nb_iter = 0                                                        #SOLUTION
        else:                                                                  #SOLUTION
            nb_iter += 1                                                       #SOLUTION
                                                                               #SOLUTION
        # on passe au meilleur voisin non tabou trouvé                         #SOLUTION
        solution_courante = meilleure                                          #SOLUTION
                                                                               #SOLUTION
        # on met à jour la liste tabou                                         #SOLUTION
        liste_tabou.append(solution_courante)                                  #SOLUTION
                                                                               #SOLUTION
    return meilleure_globale                                                   #SOLUTION

Relançons notre méthode !

nb_objets = 100
capacite = 20
random.seed(a=3)
poids_objets, valeur_objets = random_objets(10, 10)
sac = (False,)*nb_objets

print("tabou de taille 5")
sol = recherche_tabou(sac, taille_tabou=5, iter_max=30)
print("valeur finale = " + str(valeur_contenu(sol)) + ", capacite="+str(poids_contenu(sol)) + "/" + str(capacite))
print([i for i, val in enumerate(sol) if val]) # composition de la solution

tabou de taille 5
valeur finale = 62, capacite=20/20
[24, 27, 30, 41, 44, 45, 62]

2. Comportement de l'algorithme

Avec notre algorithme, on obtient une solution réalisable, mais on ne sait pas vraiment comment l'algorithme l'a obtenue. Difficile dans ces conditions d'imaginer ce qu'on pourrait faire pour améliorer la solution obtenue. Essayons déjà de voir comment l'algorithme arrive à ce résultat.

2.1 Étude de la trajectoire

Commençons par étudier la trajectoire que suit notre algorithme. C'est assez simple, il suffit de retenir toutes les solutions courantes et l'historique des meilleures solutions trouvées à chaque itération, et d'afficher tout ça avec matplotlib.

Explications Si l'on modifie le code de l'algorithme pour retenir toutes les solutions courantes et l'historique des meilleures solutions trouvées à chaque itération, et les renvoyer avec la meilleure solution sous la forme d'un tuple, on peut dessiner la trajectoire suivie par l'algorithme, et notamment toutes les dégradations qu'il a consenties au long de son exécution (cet affichage n'aura de sens que pour les méthodes par trajectoire, il faudrait trouver autre chose à afficher si on avait implémenté par exemple une méthode par population).

On va donc avoir deux séquences de données : une pour l'historique de la solution courante, et une pour l'historique de la meilleure solution trouvée depuis le début.

import matplotlib.pyplot as plt

def recherche_tabou_traject(solution_initiale, taille_tabou, iter_max):
    """
    1. On part d'un élément de notre ensemble de recherche qu'on déclare élément courant
    2. On considère le voisinage de l'élément courant, on choisit le meilleur d'entre
       eux comme nouvel élément courant parmi ceux absents de la liste tabou, et on l'ajoute
       à la liste tabou
    3. On boucle jusqu'à condition de sortie.
    """
    nb_iter = 0
    liste_tabou = deque((), maxlen = taille_tabou)

    # variables solutions pour la recherche du voisin optimal non tabou
    solution_courante = solution_initiale
    meilleure = solution_initiale
    meilleure_globale = solution_initiale

    # variables valeurs  pour la recherche du voisin optimal non tabou
    valeur_meilleure = valeur_contenu(solution_initiale)
    valeur_meilleure_globale = valeur_meilleure

    # variables pour l'affichage
    meilleure_trouvee = 0
    meilleure_globale_trouvee = 0
    
    # liste des solutions courantes et des meilleures trouvées, pour afficher la trajectoire
    # l'élément à la ième position correspond à l'itération i
    courantes =\
        deque(()) #SOLUTION
    meilleures_courantes =\
        deque(()) #SOLUTION
    
    while (nb_iter < iter_max):
        nb_iter += 1
        valeur_meilleure = 0
        
        # on parcourt tous les voisins de la solution courante
        for voisin in voisinage(solution_courante):
            valeur_voisin = valeur_contenu(voisin)                             #SOLUTION
                                                                               #SOLUTION
            # meilleure solution non taboue trouvée                            #SOLUTION
            if valeur_voisin > valeur_meilleure and voisin not in liste_tabou: #SOLUTION
                meilleure_trouvee += 1                                         #SOLUTION
                valeur_meilleure = valeur_voisin                               #SOLUTION
                meilleure = voisin                                             #SOLUTION
        
        # on met à jour la meilleure solution rencontrée depuis le début
        if valeur_meilleure > valeur_meilleure_globale:
            meilleure_globale_trouvee += 1
            meilleure_globale = meilleure
            valeur_meilleure_globale = valeur_meilleure
            nb_iter = 0
        
        meilleures_courantes.append(
            valeur_meilleure_globale) #SOLUTION
        courantes.append(
                valeur_meilleure)     #SOLUTION
        
        # on passe au meilleur voisin non tabou trouvé
        solution_courante = meilleure
        
        # on met à jour la liste tabou
        liste_tabou.append(solution_courante)   #SOLUTION

    return meilleure_globale, courantes, meilleures_courantes

sac = (False,)*nb_objets
    
iter_max = 30
taille_tabou = 5
sol, courants, meilleurs_courants=recherche_tabou_traject(sac, taille_tabou, iter_max)

plt.xlabel("nb itérations", fontsize=16)
plt.ylabel("valeur", fontsize=16)
res = plt.plot(range(len(courants)), courants)
res = plt.plot(range(len(courants)), meilleurs_courants)

On observe un phénomène intéressant au bout d'un certain temps. Comment s'appelle-t-il ? Comment l'interpréter, et que pourrions-nous faire face à ça ?

Solution On observe des oscillations entre des solutions proches, dont l'une est meilleure que l'autre. Cette solution est un optimum local, et l'algorithme semble être prisonnier de la "vallée" de cet optimum local, et qui s'appelle son bassin d'attraction. C'est ce phénomène qui justifie d'effectuer des étapes de diversification, afin de s'échapper de ce bassin d'attraction. La conclusion est qu'il faut renforcer la diversification. La première idée qui vient à l'esprit pour remédier à ce problème consiste à agrandir augmenter la taille de la liste tabou.

Si les étudiants préfèrent augmenter le nombre d'itérations (sans toucher la taille de la liste), la trajectoire ne sera pas modifiée (si ce n'est que les oscillations dureront plus longtemps). On peut les laisser tester, et constater que cette approche ne change rien à la qualité finale de la meilleure solution.

Testons cette solution (avec une valeur de 200, par exemple) :

sac=(False,)*nb_objets
taille_tabou = 200 #SOLUTION
sol, courants, meilleurs_courants=recherche_tabou_traject(sac, taille_tabou, iter_max)

plt.plot(range(len(courants)), courants)
plt.plot(range(len(courants)), meilleurs_courants)
plt.xlabel("nb itérations", fontsize=16)
plt.ylabel("valeur", fontsize=16)

Text(0, 0.5, 'valeur')

Comment interprétez-vous cette courbe ?

Solution On arrive à échapper temporairement aux oscillations, la solution est un peu améliorée. Par contre, on remarque que lorsqu'on atteint la meilleure solution trouvée, autour de la 30^ème itération, les solutions régressent pendant environ 70 itérations, avant de revenir à la même valeur. Cela suggère qu'on finit par revenir à la même solution qui est un optimum local. On en conclut qu'on est sorti du bassin d'attraction précédent, mais qu'on est tombé dans un autre bassin, plus grand, et dont on n'arrive de nouveau plus à s'échapper.

Essayons donc d'augmenter encore la valeur.

sac=(False,)*nb_objets
taille_tabou = 400 #SOLUTION
sol, courants, meilleurs_courants=recherche_tabou_traject(sac, taille_tabou, iter_max)

plt.plot(range(len(courants)), courants)
plt.plot(range(len(courants)), meilleurs_courants)
plt.xlabel("nb itérations", fontsize=16)
plt.ylabel("valeur", fontsize=16)

Text(0, 0.5, 'valeur')

Qu'est-ce que ça donne ? Comment y remédier ?

Solution Le résultat est identique. Comme on n'a pas augmenté le nombre d'itérations, la liste tabou n'a pas le temps de se remplir, et on s'arrête prématurément. Il faut aussi augmenter le nombre d'itérations.

Essayons avec la même valeur que pour la taille de la liste.

sac=(False,)*nb_objets
iter_max     = 400 #SOLUTION
taille_tabou = 400 #SOLUTION
sol, courants, meilleurs_courants=recherche_tabou_traject(sac, taille_tabou, iter_max)

plt.plot(range(len(courants)), courants)
plt.plot(range(len(courants)), meilleurs_courants)
plt.xlabel("nb itérations", fontsize=16)
plt.ylabel("valeur", fontsize=16)

Text(0, 0.5, 'valeur')

Qu'observe-t-on ?

Solution Cette fois, la recherche itère plus de fois. Mais malgré ce grand nombre d'itérations, on n'a pas augmenté la solution.

Il faut donc pousser plus loin la diversification, mais il est probable qu'augmenter encore la taille de la liste ou le nombre d'itérations ne permette d'améliorer significativement la solution que si on met des valeurs nettement plus grandes. C'est suggéré par le fait que la solution se dégrade de plus en plus, mais finit tout de même par augmenter.

Cela risque de poser de gros problèmes, tant en termes de temps de calcul que pour la mémoire. Il faudra donc utiliser d'autres techniques que simplement augmenter la taille de la liste tabou. Cette partie sera discutée plus bas.

2.2 Qualité de la solution

Mais finalement, on ne sait pas si cette solution que notre algorithme génère est vraiment bonne. On sait qu'on l'a améliorée depuis le Hill-Climbing, et qu'on peut faire mieux, puisque le multi-start y arrive. Or, la question risque d'être difficile à trancher. Pourquoi ?

Solution La question à laquelle on souhaite répondre est : est-on loin de l'optimal ? Si on en est très proche, ça ne vaut pas le coup de se lancer dans des optimisations complexes, qui rallongeront le temps de calcul, alors que les opportunités d'améliorer la solution sont minimes.

Pour comparer avec l'optimal, il faudrait être capable de déterminer cet optimal. Or le problème est $\mathcal{NP}$-Difficile. On peut rechercher manuellement une solution optimale à l'instance générée, mais c'est un peu hasardeux, il faudrait démontrer mathématiquement l'optimalité de la solution qu'on a trouvée, et le travail serait à refaire pour chaque instance. En revanche, on peut déterminer cet optimal pour des petites instances en utilisant une méthode exacte (méthode arborescente, programmation dynamique), ou en exploitant des bibliothèques d'instances déjà résolues, sur lesquelles on pourra tester notre métaheuristique.

Mais on pourrait quand même essayer de trouver un point de comparaison, moins précis, mais qui nous apporte quand même quelques informations. Et la solution est en fait très simple, elle se retrouve dans le Workshop sur la programmation linéaire. Souvenez-vous, on avait déjà travaillé sur le problème du sac à dos. Est-ce qu'on avait réussi à le résoudre ? Pourquoi ?

Solution Non, car la résolution d'un programme linéaire par l'algorithme du Simplexe ne fonctionne qu'avec des variables de décision réelles. Or nos variables de décision sont binaires, donc entières. On a affaire à un programme linéaire en nombres entiers.

Et que produisait la méthode de résolution qu'on avait essayé de mettre en place ?

Solution Une solution non réalisable, dont la valeur était supérieure à la solution optimale du problème. En effet, le simplexe génèrera une solution ne tenant pas compte des contraintes d'intégrité des variables de décision, donc une solution moins contrainte (puisqu'une contrainte réduit l'espace des solutions admissibles), meilleure en valeur, mais non réalisable.

Le voilà notre indicateur ! Grace à ça, on va pouvoir mettre en place un plan d'expériences (inutile par contre d'appliquer directement les principes mathématiques mentionnés dans ce lien, ils sont juste là pour vous donner une idée du principe général et sont des outils rarement utilisés tel quel dans le domaine de l'évaluation d'algorithmes).

Explications En fait, à part l'étude de la trajectoire, ce Workshop lui-même est organisé comme un plan d'expériences. Par exemple, le fait de considérer la recherche simple avant de pouvoir se pencher sur le multi-start, ou la définition d'un indicateur de qualité en premier lieu.

En Recherche Opérationnelle, l'indicateur de qualité d'une solution est la plupart du temps une borne, qu'elle soit obtenue par relaxation de contraintes sur un modèle linéaire, ou de manière algorithmique. C'est le cas ici, et sans le savoir, il y a deux séances, les élèves ont fait de la relaxation continue, ce qui leur a fourni une borne supérieure (puisqu'on traite un problème de maximisation).

Une fois qu'on a obtenu une borne supérieure de l'optimal pour notre instance, si la solution qu'on a généré est éloignée de la borne, que peut-on en conclure ?

Solution Rien de particulier. La solution peut être mauvaise, mais ça peut être aussi la borne qui est éloignée de l'optimal, ou une combinaison des deux.

Et si notre solution est proche de la borne ?

Solution Dans ce cas, on peut conclure que la solution est de bonne qualité. C'est donc effectivement moins révélateur que si on pouvait comparer à l'optimal, mais ça nous donne tout de même des informations. Et si on avait d'autres bornes, on pourrait essayer d'affiner ce résultat. Parfois, les problèmes $\mathcal{NP}$-Difficiles admettent des bornes algorithmiques, c'est-à-dire des bornes qui se calculent algorithmiquement. Pour le sac à dos, ce n'est pas évident, mais pour certains problèmes comme le Voyageur de Commerce, il en existe qui sont assez faciles à calculer.

Si on considère, comme ici, un problème de maximisation, il suffit alors de prendre la plus petite borne supérieure de l'instance. Et inversement, la plus grande borne inférieure si on a un problème de minimisation.

Il ne nous reste plus qu'à reprendre le code de la programmation linéaire. Testons ça :

Explications Les structures de données qu'on utilise pour notre métaheuristique sont déjà adaptées, il n'y a quasiment rien à faire. En plus, il est possible de passer par des dictionnaires et d'utiliser `lpSum()`pour générer les sommes à intégrer au modèle du problème, ça ira bien plus vite !! Cette solution sera utile si les élèves ont de gros programmes à résoudre (ils pourront même charger des fichiers Excel ou JSON pour le faire).

# on vous donne la solution toute faite !

from pulp import *

def borne_superieure():
    objets = range(nb_objets)

    # variables
    x = LpVariable.dicts('objet', objets, 0, 1)
    
    # probleme
    prob = LpProblem("knapsack", LpMaximize)

    # fonction objective
    cost = lpSum([valeur_objets[i]*x[i] for i in objets])
    prob += cost

    # contrainte
    prob += lpSum([poids_objets[i]*x[i] for i in objets]) <= capacite

    prob.solve()
    return value(prob.objective) if (LpStatus[prob.status] == "Optimal") else None

borne = borne_superieure()
if borne is not None:
    print("borne supérieure : ", borne)
print("valeur de la solution :", str(valeur_contenu(sol)))

borne supérieure :  91.33333329999999
valeur de la solution : 68

Ce n'est pas mal, mais on peut surement faire mieux ! Typiquement en faisant du multistart. D'ailleurs, voyons ce que ça donne :

import random
random.seed(a=5)

sol_max = None
val_max = 0

# multi-start de 500 itérations
for _ in range (500):
    sac = random_solution()
    sol_courante = recherche_tabou(sac, taille_tabou, iter_max)
    val_courante = valeur_contenu(sol_courante)
    if (val_courante > val_max):
        val_max = val_courante
        sol_max = sol_courante

print("valeur finale = " + str(valeur_contenu(sol_max)))

valeur finale = 91

Comment interprétez-vous ce résultat ?

Solution Puisque les valeurs des objets sont entières, la solution construite par multistart est optimale. Vue la taille réduite du problème et le nombre important de starts, l'algorithme a de très grandes chances de converger vers la solution optimale. On verra dans la suite du Workshop que les choses vont se compliquer sur des instances plus grandes.

Souvenez-vous, au Workshop précédent, deux exécutions du multistart donnaient la même valeur de solution. La voila, l'explication !

Mais on a décidé de s'en tenir aux tabou simple. Revenons-y, et essayons d'avoir une idée plus précise de ses performances.

3. Plan d'expérience et optimisation de l'algorithme

On a réussi à construire une solution aux 2/3 de la borne (voire plus si on utilise le multi-start) de notre instance aléatoire. Mais est-ce qu'on peut généraliser ce résultat ?

Solution Non, une seule mesure n'est pas représentative. On peut avoir rencontré une instance particulièrement facile. Il va falloir étudier le comportement de cet algorithme de manière un peu plus statistique.

3.1 Qualité statistique des solutions

Pour se faire une idée plus exacte des performances de l'algorithme, il faudrait tester un grand nombre d'instances. Mais cette approche pose un gros problème vis-à-vis de l'optimisation qu'on a mise en place la dernière fois pour accélérer les calculs. Quel est ce problème ?

Solution Il faut calculaer à chaque fois la distance à la borne inférieure de l'instance considérée. On évaluera ensuite de manière statistique ces distances (le choix exact des statistiques à calculer est abordé plus loin).

Le problème du point de vue du code, c'est qu'on a mis en place une stratégie de mise en cache (mémoïsation) des fonctions de calcul de la valeur et du poids, ainsi que de la génération du voisinage d'une solution. Or, la mémoïsation de ces fonctions ne fonctionne que si l'instance ne change pas, puisqu'une instance est définie par les variables globales poids_objets, poids_objets, nb_objets, et capacite. Si par exemple on change la valeur des objets, la fonction valeur_contenu renverra le calcul effectué avant le changement de valeur, et le résultat sera faux.

Comment peut-on y remédier sans pour autant renoncer à cette optimisation, qui risque d'être utile si on doit faire un grand nombre d'essais ?

Solution À moins que les étudiants aient regardé la documentation de lru_cache, il y a de fortes chances qu'ils proposent de réimplémenter les fonctions en passant à chaque fois l'instance en paramètre. Cette approche pose plusieurs problèmes  :

• La performance de la mémoïsation baisse au fur et à mesure qu'on ajoute des paramètres d'entrée, puisque ces paramètres doivent être hashés pour pouvoir les utiliser comme clé d'accès au cache.
• Ici, ça ne fonctionnera pas en l'état, car, comme dit au dessus, les arguments qu'elle prend en paramètres sont hashés pour accéder au cache. Or, les dictionnaires, qui servent pour la valeur et le poids des objets, ne sont pas hashables. En implémentant cette approche, on obtiendrait une erreur à l'exécution.

Et il existe une solution bien plus simple : effacer le cache de la fonction. Pour cela, il y a deux solutions :

• On peut vider le cache d'une fonction mémoïsée d'un simple appel à la méthode cache_clear() de la fonction (puisqu'en Python, une fonction est un objet). Pour vider le cache de la fonctionpoids_contenu, par exemple, on exécute poids_contenu.cache_clear(). Évidemment, ça doit être fait sur toutes les fonctions mémoïsées. Le plus simple est d'implémenter cette réinitialisation du cache directement dans la fonction de génération d'instance aléatoire.
• L'autre solution, plus élégante mais plus subtile, consiste à implémenter une fonction principale recherche_tabou, et à définir toutes les autres fonctions (dont celles mémoïsées), comme des fonctions imbriquées dans cette fonction principale :

  
  def recherche_tabou(nb_objets, capacite, poids_objets, valeur_objets, #instance
                      solution_initiale, taille_tabou, iter_max): #paramètres de l'algo
      @lru_cache(maxsize=None)
      def poids_contenu(sac):
          return sum(poids_objets[i] for i, val in enumerate(sac) if val)
  
      @lru_cache(maxsize=None)
      def valeur_contenu(sac):
          …
  

  Non seulement les fonctions imbriquées sont redéclarées à chaque appel du fait de la portée de ces sous-fonctions (ce qui réinitialise leur cache), mais cela permet de convertir en paramètres de la fonction principale les variables globales définissant une instance, ce qui est plus propre en terme d'organisation du code. Et les fonctions imbriquées n'ont pas besoin d'être modifiées puisqu'elles accèdent à toutes les variables visibles depuis cette fonction principale (plus précisément, à une copie locale de ces variables,si on voulait modifier ces variables dans une fonction interne il faudrait utiliser le mot-clé nonlocal). De plus, sur de grosses instances, le fait de ne pas avoir à faire la résolution de variables globales permet de gagner un peu de temps de calcul par rapport à la solution précédente.

Par soucis de lisibilité, c'est la première solution qui est présentée ici, elle est acceptable. Mais si les étudiants sont suffisamment à l'aise avec Python, la seconde est recommandée.

Corrigeons ça tout de suite !

def random_objets(poids_max, val_max):
    """
    Cette fonction genere des objets de poids et de valeur
    aleatoires (borne par les valeurs passees en parametre).
    
    Renvoie un tuple de 2 dictionnaires (poids,valeur)
    """
    
    poids_contenu.cache_clear()  #SOLUTION
    valeur_contenu.cache_clear() #SOLUTION
    voisinage.cache_clear()      #SOLUTION

    poids_objets  = {i: random.randint(1, poids_max) for i in range(nb_objets)}
    valeur_objets = {i: random.randint(1, val_max)   for i in range(nb_objets)}
    
    return poids_objets, valeur_objets

Problème résolu ! On peut attaquer l'implémentation de notre test statistique. Reprenons la taille de tabou et le nombre d'itérations qui ont donné de bons résultats tout à l'heure. Par ailleurs, pour essayer d'augmenter la représentativité de nos résultats, on va fournir à la recherche une solution initiale aléatoire à chaque test.

Et vu qu'on recommence à avoir des calculs qui prennent du temps, on va réutiliser la barre de progression du premier Workshop. Pour cette fois, le code vous est rappelé.

from ipywidgets import IntProgress
from IPython.display import display

random.seed(a=3)
bornes = deque(())

iter_max     = 200 #SOLUTION
taille_tabou = 100 #SOLUTION

# paramètres du test
nb_tests = 200

# on affiche la barre de progression
bar = IntProgress(min=0, max=nb_tests, layout={"width" : "100%"})
display(bar)

for _ in range(nb_tests):
    # génération aléatoire de l'instance et calcul de   #SOLUTION
    # sa borne supérieure.                              #SOLUTION
    poids_contenu.cache_clear()                         #SOLUTION
    valeur_contenu.cache_clear()                        #SOLUTION
    poids_objets, valeur_objets = random_objets(10, 10) #SOLUTION
    borne = borne_superieure()                          #SOLUTION
                                                        #SOLUTION
    sac = random_solution()                             #SOLUTION
    sol = recherche_tabou(sac, taille_tabou, iter_max)  #SOLUTION
    val = valeur_contenu(sol)                           #SOLUTION
    distance = borne-val                                #SOLUTION
    bornes.append(distance)
    bar.value += 1

# on cache la barre de progression
bar.close()

print(sum(bornes) / len(bornes))

IntProgress(value=0, layout=Layout(width='100%'), max=200)

15.4633333339

Regardons un peu plus en détail les valeurs qu'on obtient. Ça tombe bien, on a la liste de toutes les distances. Comment pourrait-on examiner plus en détail ces valeurs ?

Solution On peut afficher la distribution des distances sous forme d'histogramme.

Essayons ça !

plt.hist(bornes, edgecolor = "black")                         #SOLUTION
plt.xlabel("distance à la borne")                             #SOLUTION
plt.ylabel("nombre d'intances")                               #SOLUTION
plt.title("Distribution des distances à la borne supérieure") #SOLUTION

Text(0.5, 1.0, 'Distribution des distances à la borne supérieure')

Le gros des valeurs se situe entre 10 et 30. Seulement, ça ne nous dit pas grand chose. Quelle information nous manque-t-il ?

Solution Il faudrait mettre en relation la distance avec la valeur de la borne, sous forme par exemple d'un pourcentage $\frac{\text{valeur solution}}{\text{valeur borne}}\times 100$.

Pas de chance, il faut refaire les calculs !

Explications Ici, les étudiants n'ont perdu que quelques instants, mais on peut leur faire remarquer que sur des calculs plus longs (et ils vont arriver très vite), ne pas réfléchir aux données qu'on veut obtenir avant de lancer les calculs peut leur faire perdre beaucoup de temps. Ce point est très important quand on conçoit un plan d'expérience, il faut bien réfléchir aux indicateurs les plus pertinents.

random.seed(a=3)
bornes = deque(())

# on affiche la barre de progression
bar = IntProgress(min=0, max=nb_tests, layout={"width" : "100%"})
display(bar)

for _ in range(nb_tests):
    # génération aléatoire de l'instance et calcul de   #SOLUTION
    # sa borne supérieure.                              #SOLUTION
    poids_objets, valeur_objets = random_objets(10, 10) #SOLUTION
    borne = borne_superieure()                          #SOLUTION
                                                        #SOLUTION
    sac = random_solution()                             #SOLUTION
    sol = recherche_tabou(sac, taille_tabou, iter_max)  #SOLUTION
    val = valeur_contenu(sol)                           #SOLUTION
    distance = val*100/borne                            #SOLUTION
    bornes.append(distance)                             #SOLUTION
    bar.value += 1                                      #SOLUTION

# on cache la barre de progression
bar.close()

plt.hist(bornes, edgecolor = "black")                         #SOLUTION
plt.xlabel("distance à la borne (pourcentage)")               #SOLUTION
plt.ylabel("nombre d'intances")                               #SOLUTION
plt.title("Distribution des distances à la borne supérieure") #SOLUTION

IntProgress(value=0, layout=Layout(width='100%'), max=200)

Text(0.5, 1.0, 'Distribution des distances à la borne supérieure')

La représentation est un peu grossière, matplotlib a utilisé un intervalle de regroupement des valeurs trop grand. Quel paramètre va-t-on devoir modifier pour faire un affichage plus fin ?

Solution On va jouer sur le paramètre bins pour avoir une vision un peu plus précise.

Coup de chance, cette fois-ci, pas besoin de recalculer, on a déjà les valeurs dans la variable bornes. Il suffit de refaire un affichage de ces valeurs :

plt.hist(bornes, bins=50, edgecolor = "black"); #SOLUTION

Là, c'est plus lisible ! On voit bien que la grande majorité de nos solutions se situent au dessus de 70% de la borne. Quelle est la valeur la plus représentative du comportement de l'algorithme ? Quelles statistiques descriptives pourraient répondre à cette question ?

Solution On peut regarder la valeur médiane, ou la moyenne. La moyenne semble plus intéressante, puisqu'on peut lui associer l'écart-type (même s'il est aussi possible d'associer à la médiane l'écart absolu médian, mais les étudiants ne connaissent probablement pas cette statistique), afin de se faire une idée de la dispersion des valeurs autour de la moyenne.

Regardons ce que ça donne. Selon les performances, on pourra utiliser les fonction du module de statistiques standard de Python, ou les fonctions statistiques de NumPy.

Explications Un test rapide montre les choix d'implémentation optimaux :

• Pour la moyenne, la fonction mean du module statistics a un gros défaut : elle est très lente (du fait de sa stabilité numérique, pas très utile ici). Mieux vaut se tourner vers fmean du même module, moins précise mais nettement plus rapide (et aussi plus rapide que la version de NumpPy ou qu'une implémentation manuelle du type sum(bornes)/len(bornes).
• Pour l'écart-type, la fonction std de NumPy est la plus rapide.
• Pour la médiane, c'est de nouveau le module statistics qui propose la version la plus rapide avec la fonction median.

Il est possible de fournir directement ces résultats aux étudiants.

import statistics
import numpy as np

print("médiane={:.2f}".format(statistics.median(bornes))) #SOLUTION
print("moyenne={:.2f}, écart-type={:.2f}".format(         #SOLUTION
      statistics.fmean(bornes),np.std(bornes)))           #SOLUTION

médiane=84.76
moyenne=84.37, écart-type=7.89

Comment interprétez-vous ces résultats ?

Solution Ces valeurs sont assez satisfaisantes. Le fait que la moyenne et la médiane soient proches s'explique par le fait que la distribution est assez symétrique. On va donc garder le couple moyenne/écart-type.

On a donc un comportement plutôt satisfaisant de l'algorithme pour le type d'instances qu'on considère.

3.2 Optimisation des paramètres de l'algorithme

Si le comportement semble assez acceptable, il y a tout de même des instances pour lesquelles les choses se passent un peu moins bien. Est-ce qu'on ne pourrait pas trouver un paramétrage de l'algorithme qui s'en sort mieux avec ces instances ?

Explications On a bien vu tout à l'heure que changer la taille de la liste tabou jouait, mais qu'il fallait ajuster le nombre d'itération en fonction. La question qu'on peut se poser, c'est : quel nombre d'itération est optimal pour une taille donnée de la liste tabou ? Et maintenant qu'on sait mesure la qualité d'une solution, on va pouvoir répondre à cette question.

Pour répondre à cette question, on va regarder l'impact de la taille de la liste tabou sur la qualité. Le but va être d'afficher avec Matplotlib la moyenne sous la forme d'une courbe, et l'écart-type par une "bande" entourant la courbe :

Explications On va fixer un nombre d'itérations, faire varier la taille de la liste tabou, et calculer pour chaque configuration (chaque combinaison de nombre d'itérations et de taille de liste) la qualité moyenne et son écart-type.

Pour afficher la bande correspondant à l'écart-type, on va utiliser la fonction fill_between (déjà utilisée pour la résolution graphique d'un programme linéaire). Les limites de la surface à afficher seront simplement obtenues en additionnant (limite haute) ou soustrayant (limite basse) les écarts-types aux valeurs moyennes. Le code pour l'affichage est fourni.

Pour le reste, le principe est quasiment le même qu'avant, il faudra juste faire une double boucle imbriquée. On repart du premier paramétrage du nombre d'itérations (20) pour aller un peu plus vite, et on fait varier la taille de la liste tabou de 1 à 20. Les calculs étant un peu plus longs, on va baisser le nombre de tests par configuration à 100.

Par contre, petit détail : comme cette fois on va avoir une double boucle imbriquée (d'où le temps de calcul rallongé), il va falloir calculer le nombre d'étapes de la barre de progression !

random.seed(a=3)

# paramètres du test
tabou_min = 1
tabou_max = 20
nb_objets = 200

nb_test  = 100
iter_max = 20

# on affiche la barre de progression
nb_steps_bar =\
    (tabou_max-tabou_min)*nb_test #SOLUTION
bar = IntProgress(min=0, max=nb_steps_bar, layout={"width" : "100%"})
display(bar)

# pour stocker les résultats
moyennes   = []
deviations = []

random.seed(a=3)

# cette fois on boucle sur la taille de la liste tabou
for taille_tabou in range(tabou_min, tabou_max):            #SOLUTION
    bornes = deque(())                                      #SOLUTION
    for _ in range(nb_test):                                #SOLUTION
        # génération aléatoire de l'instance et calcul de   #SOLUTION
        # sa borne supérieure.                              #SOLUTION
        poids_objets, valeur_objets = random_objets(10, 10) #SOLUTION
        borne = borne_superieure()                          #SOLUTION
                                                            #SOLUTION
        sac = random_solution()                             #SOLUTION
        sol = recherche_tabou(sac, taille_tabou, iter_max)  #SOLUTION
        val = valeur_contenu(sol)                           #SOLUTION
        bar.value += 1                                      #SOLUTION
        distance = val*100/borne                            #SOLUTION
        bornes.append(distance)                             #SOLUTION
                                                            #SOLUTION
    # en fait, pour la moyenne, le calcul manuel reste      #SOLUTION
    # l'implémentation la plus rapide                       #SOLUTION
    moyennes.append(statistics.fmean(bornes))               #SOLUTION
    deviations.append(np.std(bornes))                       #SOLUTION

# on cache la barre de progression
bar.close()

# affichage de la courbe de moyenne
plt.plot(range(tabou_min, tabou_max), moyennes)

# affichage de la bande d'écart-type
plt.fill_between(range(tabou_min, tabou_max),
                 np.subtract(moyennes, deviations), # borne haute
                 np.add(moyennes, deviations),      # borne basse
                 alpha=.1)                          # transparence
plt.xlabel("taille de la liste tabou")
plt.ylabel("distance à la borne")
plt.title("Impact de la taille de la liste tabou sur la qualité des solutions")
plt.show()

IntProgress(value=0, layout=Layout(width='100%'), max=1900)

Qu'observez-vous ?

Solution On observe le phénomène qu'on avait déjà rencontré pus haut. Quand on augmente la taille de la liste tabou, la qualité des solutions s'améliore. Mais on finit par atteindre une situation stationnaire. Comme on l'avait déjà vu, si on ne change pas le nombre d'itérations de la recherche, à partir d'un certain niveau la liste n'a pas le temps d'être remplie. Comme ce nombre d'itérations est trop faible, agrandir la taille de la liste tabou ne change en fait rien à partir d'une certaine valeur.

Ici, la situation se stabilise vers une taille de 10. En dessous, on générera des solutions de moins bonne qualité. Au dessus, on rallongera considérablement le temps de calcul sans gain notable.

À noter, on obtient le même profil général si on augmente le nombre de tests, c'est donc un résultat représentatif.

Cela ouvre des perspectives pour continuer les tests, et en tirer des conclusions sur les paramètres les plus efficaces de l'algorithme. Comment pourrait-on procéder ?

Solution Il faudrait tester d'autres combinaisons taille de la liste tabou / nombre d'itérations de la recherche, afin de trouver les combinaisons les plus efficaces. Idéalement, on en déduirait une définition fonctionnelle optimale de ces paramètres (la valeur d'un des paramètres permet de calculer l'autre).

Voici la courbe qu'on obtient lorsqu'on revient à 200 itérations, en faisant varier la taille de la liste tabou entre 10 et 200 (quelques heures de calcul), et avec 2000 tests par itération pour avoir un résultat plus représentatif. Affichez-la et discutez avec les étudiants de ce résultat.

On remarque que dans les deux cas, une taille de liste obtenue en divisant par 2 le nombre d'itérations semble être le bon paramétrage, en tous cas pour le type d'instances qu'on a considéré. C'est évidemment à vérifier avec d'autres expériences.

Si on veut vérifier ça de manière automatique, comment pourrait-on procéder ?

Solution Dans la mesure où la qualité ne se stabilise jamais complètement, il faut décider d'un seuil de variation en dessous duquel on considère que la solution ne s'améliore plus (sous forme de pourcentage). On pourra ensuite déterminer automatiquement l'évolution de la valeur de stabilisation de la taille de la liste tabou, en fonction de l'évolution du nombre d'itérations. On pourra vérifier si notre hypothèse d'une relation linéaire de coefficient 2 entre ces deux paramètres se vérifie.

En tous cas, tout porte à croire que notre paramétrage de l'algorithme n'était pas optimal. On perdait en temps de calcul par rapport à la qualité de la solution générée.

3.3 Impact des paramètres d'instance

Mais finalement, on veut améliorer les performances sur ces instances difficiles, alors qu'on ne sait même pas ce qui les rend compliquées à résoudre. Du coup, la question qu'on peut se poser, c'est : qu'est-ce qui rend une instance difficile ?

Explications Tous les paramètres d'un instance peuvent être impactants sur le comportement de l'algorithme : * La taille du sac * Le nombre d'objets disponibles * La taille des objets (taille moyenne, mais aussi écart-type, voire distribution)

Et plus que leurs valeurs absolues, c'est leur rapport qui va jouer. Lorsqu'on aura examiné ces paramètres, on pourra continuer en étudiant :

• Le rapport entre la capacité du sac et le nombre d'objets
• Le rapport entre la capacité du sac et la taille des objets (taille moyenne, écart-type)
• Une corrélation (par exemple relation linéaire) entre le poids des objets et leur valeur. Certaines études ont montré que la présence d'une telle relation pouvait changer le comportement des algorithmes d'optimisation.

Pour commencer, on va jouer sur le nombre d'objets, c'est ce qu'il y a de plus simple à implémenter.

Explications Pour avoir des résultats représentatifs, il faudra, pour chaque nombre d'objets, calculer une qualité moyenne des solutions de l'algorithme. Par conséquent, on devra associer le nombre d'objets avec la distance moyenne à la borne (toujours en pourcentage). On pourrait utiliser une liste de tuples `(taille, distance)`, mais un dictionnaire `{taille: distance}` fera la même chose de manière plus simple et efficace. On pourra fournir les bonnes valeurs à afficher à Matplotlib avec les méthodes keys() et values().

Attention, ici on garde 200 tests pour chaque configuration (chaque valeur de nombre d'objets), mais le calcul risque de prendre du temps (même avec la mémoïsation). Plus encore qu'avant, certains élèves auront probablement besoin de réduire le nombre d'instances générées pendant le Workshop. On va aussi introduire un paramètre objets_steps, qui nous permettra de ne pas considérer toutes les valeurs de nombre d'objets. Les courbes seront plus hachées, mais les résultats n'en seront pas moins représentatifs pour autant. Si les calculs prennent trop de temps, on peut augmenter la taille des steps. Par contre, le calcul du nombre d'étapes de la barre de progression doit être adapté.

Et toujours pour gagner du temps, on reprend la valeur initiale de 5 pour la liste tabou pour gagner du temps. Grace à nos travaux précédents, on sait adapter le nombre d'itération itérations à cette taille de liste.

random.seed(a=3)

# paramètres du test
objets_min   = 20
objets_max   = 220
objets_steps = 20
nb_test      = 200

taille_tabou = 5
iter_max     = taille_tabou*2

nb_steps_bar =\
    ((objets_max-objets_min)/objets_steps)*nb_test #SOLUTION

# on affiche la barre de progression
bar = IntProgress(min=0, max=nb_steps_bar, layout={"width" : "100%"})
display(bar)

# pour stocker les résultats
moyennes   = []
deviations = []

# on fait varier le nombre d'objets dans la génération aléatoire d'instances 
for nb_objets in range(objets_min, objets_max, objets_steps): #SOLUTION
    bornes = deque(())                                        #SOLUTION
    for _ in range(nb_test):                                  #SOLUTION
        # génération aléatoire de l'instance et calcul de sa  #SOLUTION
        # borne supérieure.                                   #SOLUTION
        poids_objets, valeur_objets = random_objets(10, 10)   #SOLUTION
        borne = borne_superieure()                            #SOLUTION
                                                              #SOLUTION
        sac = random_solution()                               #SOLUTION
        sol = recherche_tabou(sac, taille_tabou, iter_max)    #SOLUTION
        val = valeur_contenu(sol)                             #SOLUTION
        bar.value += 1                                        #SOLUTION
        distance = val*100/borne                              #SOLUTION
        bornes.append(distance)                               #SOLUTION
    moyennes.append(statistics.fmean(bornes))                 #SOLUTION
    deviations.append(np.std(bornes))                         #SOLUTION

# on cache la barre de progression
bar.close()

# affichage de la courbe de moyenne
plt.plot(range(objets_min, objets_max, objets_steps), moyennes)

# affichage de la bande d'écart-type
plt.fill_between(range(objets_min, objets_max, objets_steps),
                 np.subtract(moyennes, deviations), # borne haute
                 np.add(moyennes, deviations),      # borne basse
                 alpha=.1)                          # transparence
plt.ylim(top=100)
plt.xlabel("nombre d'objets disponibles")
plt.ylabel("distance à la borne")
plt.title("Impact du nombre d'objets sur la qualité des solutions")
plt.show()

IntProgress(value=0, layout=Layout(width='100%'), max=2000)

Que pensez-vous de ces résultats ? Comment les interprétez-vous ?

Ces résultats sont tout à fait prévisibles. Plus le problème est grand, plus le nombre de combinaisons d'objets est important, ce qui rend la recherche d'une bonne solution plus complexe (ce qui explique aussi que l'écart-type augmente). À première vue, la baisse de qualité semble être de moins en moins forte au fur et à mesure qu'on augmente le nombre d'objets disponibles. On peut se demander si ça continuerait à baisser à la même cadence si on augmentait encore le nombre d'objets. Pour cela, il faut augmenter la taille maximale du nombre d'objets, quitte à prendre des steps plus grands pour ne pas perdre trop de temps.

Voici les courbes qu'on obtient si on teste cette approche, et en faisant 500 tests par configuration, pour les deux configurations de l'algorithme qu'on a utilisées jusqu'ici, obtenue au prix d'un calcul de plusieurs heures (plusieurs jours pour le paramétrage à 200 itérations). Affichez-la et commentez-la avec les étudiants.

Les deux courbes ont quasiment le même profil, mais l'algorithme se comporte mieux lorsqu'on augmente la taille de la liste tabou et le nombre d'itérations (attention à la différence d'échelle entre les deux graphiques). Et lorsque le nombre d'objets augmente, la qualité diminue de plus en plus lentement. Dans le cas du paramètre à 10 itérations, la qualité finit par se stabiliser. En revanche, dans le cas de 200 itération, la qualité remonte légèrement.

Une des interprétations possibles de ces résultats est qu'en augmentant le nombre de variables, on offre à l'algorithme plus d'opportunités de s'échapper d'un puits d'attraction (ce qui pourrait même expliquer la légère remontée dans le cas de 200 itérations).

Il ne faut pas non plus perdre de vue que la distance ne se calcule pas à l'optimal mais à une borne. Il est possible que la qualité de cette borne évolue au cours du temps, et fausse ces résultats.

Et puis, si on veut aller encore plus loin, on peut essayer d'identifier les types d'instances qui posent problèmes à l'algorithme (les instances dites "difficiles"), et examiner spécifiquement ces instances. Elles ne sont pas compliquées à retrouver, ce sont celles qui apparaissent sur les barres les plus à gauche dans vos histogrammes de distance à la borne. Si ça se trouve, elles ont des caractéristiques particulières en terme de nombre, valeur ou de poids d'objets, de taille de sac…

L'idéal, ce serait d'être capable de fixer les paramètres de l'algorithme (taille de liste tabou, nombre d'itérations) en fonction des paramètres de l'instance. Cette approche fonctionnelle des paramètres de l'algorithme s'appelle optimisation online, par opposition à l'optimisation offline qui les fixe de manière absolue (indépendant de la forme de l'instance). On peut même utiliser une méthode d'optimisation pour déterminer les paramètres les plus efficaces, en général en introduisant un aspect aléatoire. Cette technique s'appelle métaoptimisation.

3.4 Pour aller plus loin

Il y a encore beaucoup de choses qu'on pourrait étudier pour comprendre et surtout améliorer le fonctionnement de notre algorithme (sans compter les calculs de vous devriez peut-être refaire avec un plus grand nombre de tests par configuration). On n'a regardé que des paramètres seuls, pas la manière dont ils peuvent se combiner (ou très peu). Et puis on ne s'est pas du tout penchés sur le temps de calcul. Comme lorsqu'on étudiait la complexité théorique des algorithmes, mesurer le temps d'occupation CPU ne sera pas forcément pertinent. Quelle métrique vous semble avoir du sens ?

Solution Le nombre d'itérations nécessaire pour se stabiliser (en utilisant le principe de seuil évoqué plus haut) est un bon indicateur. Attention, lorsqu'on parle du multi-start, on peut compter uniquement le nombre d'itérations totale de toutes les recherches effectuées, ou compter aussi le nombre de starts.

Malgré tout, il est parfois nécessaire de compter le temps CPU, car ce nombre d'itération n'est pas forcément représentatif. On a pu voir de manière empirique, lors du test faisant varier le nombre d'objets, que les instances de petite taille étaient traitées plus rapidement, ce qui s'explique en partie par la taille du voisinage à traiter à chaque itération. Or, cette taille n'impacte pas directement le nombre d'itérations.

En tous cas, ça nous a permis d'améliorer la qualité de la solution en affinant les paramètres de la recherche tabou pour renforcer la phase de diversification (qui ici semblait nous poser problème). Si on rencontrait des instances plus difficiles, que notre algorithme n'arrive pas à résoudre efficacement, quel que soit le paramétrage, que pourrions-nous faire ?

Solution Dans ce cas, c'est le principe même de l'étape de diversification qui doit être amélioré.

Parmi les techniques existantes, on peut citer :

• On peut utiliser une liste tabou dynamique, dont le comportement est plus évolué (typiquement, on réduit la taille lorsque la solution courante s'améliore, on l'augmente lorsqu'elle diminue, et on la vide lorsque la meilleure solution globale s'améliore).
• On peut utiliser un voisinage plus adapté. Par exemple, pour favoriser la diversification, on peut considérer comme voisin toutes les solutions obtenues en inversant un ou deux objets.
• On peut agir sur la génération du voisinage, par exemple en ne générant qu'un sous-ensemble des voisins, de manière à diriger la recherche (méthode par liste de candidats).
• On peut accepter de passer par des solutions non réalisables (appelées structures de référence) pour atteindre de meilleures solutions réalisables (méthode par chaînes d'éjection).
• On peut implémenter des mécanismes de mémoire à long terme en plus de la liste tabou (mémoire à court terme). Cette approche nécessite de modéliser la notion de type de mouvement permettant de passer d'une solution à l'autre, en plus (voire à la place) de la notion de voisinage. Ici, un mouvement pourrait consister à ajouter un objet en particulier, ou à ajouter deux objets conjointement, alors qu'un voisin est constitué d'une solution complète. On favoriserait alors les mouvements ayant mené à de meilleures solutions.
• Etc. (il existe un très grand nombre de sophistications de la méthode tabou, y compris pour le problème du sac à dos, la littérature à ce sujet est très abondante).

Certaines de ces techniques sont parfois assez complexes, et nécessitent beaucoup d'expérimentation pour donner de bons résultats pour un problème donné.

Pour finir, on peut aussi choisir d'implémenter une métaheuristique différente, mais de toute façon ce travail d'optimisation du paramétrage sera nécessaire.

Et puis une fois tout ce travail fait (et pas avant), on pourrait aussi regarder le multi-start. Qu'est-ce qui vous semble pertinent ?

Le seul paramètre de l'algorithme qu'on peut regarder, c'est le nombre de starts. Évidemment, plus on en fait mieux c'est. Mais le problème, ce sera le temps de calcul. La question qui se pose, c'est est-il plus intéressant de faire des recherches optimisées mais qui prennent du temps, ou des recherches moins efficaces mais plus rapides, qui permettent donc d'effectuer plus de starts ? Pour y répondre, il faut tester des combinaisons (nombre de starts,nombre d'itérations de la recherche), en utilisant à chaque fois la taille de liste tabou la plus adaptée, en s'appuyant sur les résultats précédents.

Le cas extrême est celui d'un Hill-Climbing multi-start. Un test rapide permettra en quelques minutes de vérifier qu'il ne donnera de meilleurs résultats qu'un tabou simple qu'au prix d'un nombre très important de starts. Finalement, cela prendra plus de temps pour avoir le même résultat.

Là aussi, on pourra en fonction de ces résultats optimiser le paramétrage du multi-start. Non seulement en définissant le nombre de starts, mais aussi en améliorant la génération aléatoire de solution de départ. On peut par exemple utiliser la mémoire à long terme dont il est question au dessus pour orienter la génération de cette solution de départ.

4 Conclusion

Ce Workshop est maintenant terminé, et on a une idée beaucoup plus précise du comportement de notre métaheuristique.

Explications On a vu que le paramétrage fin des différentes étapes d'une métaheuristique peut grandement impacter son comportement, et le côté stochastique est loin d'être trivial, particulièrement lorsqu'un aspect aléatoire est présent dans l'algorithme (ce qui est par exemple le cas pour le multi-start et le GRASP, ou lorsqu'on fait de la métaoptimisation). La génération de solutions aléatoires de manière uniforme est lui-même un sujet assez compliqué, et qui fait intervenir des notions avancées de probabilités et de processus stochastiques.

En fait, ce Workshop lui-même est organisé comme un plan d'expériences (par exemple, le fait de considérer la recherche simple avant de pouvoir se pencher sur le multi-start, ou la définition d'un indicateur de qualité en premier lieu).

Il faudra réfléchir à un plan similaire pour votre projet, de manière à ne pas perdre de temps (d'autant que les calculs sont plutôt longs). Par contre, ne poussez pas la démarche jusqu'à appliquer les principes mathématiques mentionnés dans le lien précédent, ce sont des outils rarement utilisés tel quel dans le domaine de l'évaluation d'algorithmes.

Explications La notion de plan d'expériences est fondamentale, il faut vraiment sensibiliser les étudiants au fait que l'ensemble de ce Workshop est l'implémentation d'un plan d'expérience (à part une erreur volontaire sur l'indicateur de qualité, qui sert à illustrer le coût d'une erreur de conception). C'est l'un des livrables de leur projet, et un plan d'expérience mal conçu ne leur permettra pas d'obtenir des résultats révélateurs, notamment à cause du temps de calcul prohibitif.

Il faut par contre insister sur le fait que les méthodes formelles issues des plans d'expériences numériques n'est pas du tout à maîtriser ici. Par ailleurs, il est fortement recommandé de stocker les résultats des calculs (par exemple dans des fichiers JSON) avec une convention de nommage qui permet de retrouver facilement le paramétrage de l'expérience.

Vous avez maintenant toutes les armes pour, non seulement concevoir et implémenter une métaheuristique, mais aussi analyser son comportement et optimiser ses paramètres pour en tirer le meilleur.