Planifiez vos expériences avec un plan factoriel

Le script nocturne d’agrégation des ventes est trop lent. L’équipe a trois idées : lire les données au format Parquet plutôt qu’en CSV, activer un cache des tables de référence, et paralléliser le calcul sur 8 processus. Chaque essai prend une dizaine de minutes. Faut-il tester les trois idées une par une ? Toutes les combinaisons ? Et comment savoir laquelle compte vraiment ?

C’est l’objet des plans d’expériences (design of experiments) : obtenir le maximum d’information avec le minimum d’essais, sans passer à côté des interactions entre facteurs.

Trois façons de planifier

Rappelons le vocabulaire du chapitre 1 : un facteur est un paramètre que l’on fait varier (le format, le cache, le parallélisme), un niveau est une valeur de ce facteur (CSV ou Parquet), la réponse est la métrique mesurée (ici, la durée du script).

Plan	Nombre d’essais pour $k$ facteurs à $n_i$ niveaux	Avantage	Limite
Un facteur à la fois	$1 + \sum_i (n_i - 1)$	très peu d’essais	ne voit pas les interactions, peut conclure à tort
Factoriel complet	$\prod_i n_i$	toute l’information	explose avec le nombre de facteurs
Plan $2^k$	$2^k$	voit toutes les interactions, analyse simple	deux niveaux seulement par facteur
Plan fractionnaire $2^{k-p}$	$2^{k-p}$	encore moins d’essais	certains effets sont confondus

Avec 3 facteurs à 4 niveaux, le plan « un facteur à la fois » demande $1 + 3 \times 3 = 10$ essais, le factoriel complet $4^3 = 64$.

Avec trois facteurs à deux niveaux : un facteur à la fois ne s'éloigne jamais du point de départ, le plan 2³ couvre tout le cube, le plan fractionnaire en couvre la moitié, choisie pour rester équilibrée.

Pourquoi le plan « un facteur à la fois » est-il dangereux ? Partez de la configuration « CSV, sans cache, 1 processus », qui prend 120 s. Passer à 8 processus fait gagner 50 s ; passer à Parquet en fait gagner 75. Vous en déduisez que les deux ensemble feront gagner 125 s… sur un script de 120 s ! Le gain réel est bien moindre : une fois la lecture accélérée par Parquet, il reste beaucoup moins de travail à paralléliser. Deux facteurs interagissent quand l’effet de l’un dépend du niveau de l’autre, et c’est précisément ce qu’un plan « un facteur à la fois » ne peut pas voir.

Le plan $2^2$ : effets et interactions

On commence simple, avec deux facteurs à deux niveaux, codés $-1$ et $+1$ :

• $A$ : le format ($-1$ = CSV, $+1$ = Parquet) ;
• $B$ : le cache ($-1$ = désactivé, $+1$ = activé).

Les quatre essais donnent les durées suivantes (en secondes) :

Essai	$x_A$	$x_B$	Durée $y$
1	−1	−1	80
2	+1	−1	30
3	−1	+1	60
4	+1	+1	20

On cherche un modèle de la forme :

$$y = q_0 + q_A x_A + q_B x_B + q_{AB} x_A x_B$$

• $q_0$ est la réponse moyenne ;
• $q_A$ est l’effet du facteur $A$ : passer de $-1$ à $0$ (ou de $0$ à $+1$) change $y$ de $q_A$ ;
• $q_{AB}$ est l’interaction entre $A$ et $B$.

Quatre équations, quatre inconnues : le système a une solution unique. Et grâce au codage en $\pm 1$, la solution est très simple à calculer avec la table des signes.

La table des signes

On écrit une colonne par terme du modèle. La colonne $I$ ne contient que des $+1$ ; la colonne $AB$ est le produit des colonnes $A$ et $B$. Chaque effet s’obtient en multipliant la colonne $y$ par la colonne de signes correspondante, en sommant, puis en divisant par $2^k = 4$.

$I$	$A$	$B$	$AB$	$y$
+1	−1	−1	+1	80
+1	+1	−1	−1	30
+1	−1	+1	−1	60
+1	+1	+1	+1	20
190	−90	−30	10	total
47,5	−22,5	−7,5	2,5	total / 4

Donc $q_0 = 47{,}5$ s, $q_A = -22{,}5$ s, $q_B = -7{,}5$ s et $q_{AB} = 2{,}5$ s. Le modèle est $y = 47{,}5 - 22{,}5 x_A - 7{,}5 x_B + 2{,}5 x_A x_B$. Vérifions l’essai 1 : $47{,}5 + 22{,}5 + 7{,}5 + 2{,}5 = 80$. Parfait.

Passer à Parquet fait gagner 50 s sans cache et 40 s avec : les pentes diffèrent peu, d'où une interaction faible, de 2,5 s. Des droites parallèles signifieraient une interaction nulle.

Quelle part de la variation chaque facteur explique-t-il ?

L’effet brut ne suffit pas à comparer les facteurs : on veut savoir quelle part de la variation totale de $y$ chacun explique. La variation totale est $\text{SST} = \sum_j (y_j - \bar{y})^2$, et l’on montre qu’elle se décompose exactement en :

$$\text{SST} = 2^k q_A^2 + 2^k q_B^2 + 2^k q_{AB}^2 = \text{SS}_A + \text{SS}_B + \text{SS}_{AB}$$

Ici, $\text{SS}_A = 4 \times 22{,}5^2 = 2\,025$, $\text{SS}_B = 225$, $\text{SS}_{AB} = 25$, et $\text{SST} = 2\,275$. Le format explique 89 % de la variation, le cache 9,9 %, leur interaction 1,1 %. Si l’on ne peut faire qu’une chose, c’est passer à Parquet.

À retenir

Dans un plan $2^k$ codé en $\pm 1$, chaque effet vaut $q = \frac{1}{2^k} \sum_j s_j y_j$, où $s_j$ est le signe de la colonne correspondante. La part de variation expliquée par un effet est $2^k q^2 / \text{SST}$. On classe les facteurs par cette part, pas par intuition.

Le plan $2^3$, en Python

Ajoutons le troisième facteur, $C$ : le parallélisme ($-1$ = 1 processus, $+1$ = 8 processus). Il faut maintenant $2^3 = 8$ essais, et le modèle comporte 8 termes : la moyenne, 3 effets principaux, 3 interactions doubles ($AB$, $AC$, $BC$) et une interaction triple ($ABC$). Le calcul à la main devient fastidieux ; laissons Python faire :

import itertools

# Ordre standard : A varie le plus vite, puis B, puis C
essais = [(a, b, c) for c, b, a in itertools.product([-1, 1], repeat=3)]
y = [120, 45, 100, 35, 70, 30, 55, 25]          # durées mesurées (s)

colonnes = {
    "A": lambda a, b, c: a,      "B": lambda a, b, c: b,      "C": lambda a, b, c: c,
    "AB": lambda a, b, c: a * b, "AC": lambda a, b, c: a * c, "BC": lambda a, b, c: b * c,
    "ABC": lambda a, b, c: a * b * c,
}
n = len(y)
effets = {nom: sum(f(*e) * yj for e, yj in zip(essais, y)) / n for nom, f in colonnes.items()}
SST = sum(n * q ** 2 for q in effets.values())

print(f"q0 = {sum(y) / n}")
for nom, q in sorted(effets.items(), key=lambda kv: -abs(kv[1])):
    print(f"{nom:<4} q = {q:7.2f}   part de variation = {100 * n * q ** 2 / SST:5.1f} %")

Le résultat est parlant :

Effet	$q$ (s)	Part de variation
$A$ (format)	−26,25	66,4 %
$C$ (parallélisme)	−15,00	21,7 %
$AC$	+8,75	7,4 %
$B$ (cache)	−6,25	3,8 %
$AB$	+2,50	0,6 %
$BC$	+1,25	0,2 %
$ABC$	0	0 %

Le format domine, suivi du parallélisme. L’interaction $AC$ n’est pas négligeable (7,4 %) et elle est positive alors que $q_A$ et $q_C$ sont négatifs : les deux améliorations ne s’additionnent pas complètement. Le parallélisme fait gagner 50 s en CSV sans cache (de 120 à 70 s) mais seulement 15 s en Parquet sans cache (de 45 à 30 s). Une fois la lecture accélérée, il reste moins de calcul à paralléliser.

Trois termes expliquent 95 % de la variation : le format, le parallélisme et leur interaction AC. À droite, l'interaction se voit : le parallélisme fait gagner 50 s en CSV mais 15 s seulement en Parquet.

Attention

Le signe d’une interaction ne se lit pas seul. Ici, $q_{AC} > 0$ signifie que l’effet de $C$ est moins négatif quand $A = +1$, donc que le parallélisme rapporte moins une fois qu’on est passé à Parquet. Revenez toujours aux durées mesurées pour interpréter.

Répéter pour estimer l’erreur : le plan $2^k r$

Jusqu’ici, chaque configuration n’a été mesurée qu’une fois. Or le chapitre 7 l’a montré : une mesure est une variable aléatoire. Comment savoir si un effet de 2,5 s est réel ou s’il n’est que du bruit ? En répétant chaque essai $r$ fois : c’est le plan $2^k r$.

Reprenons le plan $2^2$ avec $r = 3$ répétitions :

Essai	$x_A$	$x_B$	Mesures (s)	Moyenne
1	−1	−1	78, 82, 80	80
2	+1	−1	29, 31, 30	30
3	−1	+1	57, 61, 62	60
4	+1	+1	21, 19, 20	20

Les effets se calculent sur les moyennes, comme avant : $q_A = -22{,}5$, $q_B = -7{,}5$, $q_{AB} = 2{,}5$. La nouveauté est l’erreur expérimentale : l’écart de chaque mesure à la moyenne de son essai. Sa somme des carrés vaut $\text{SSE} = 4 + 4 + 0 + 1 + 1 + 0 + 9 + 1 + 4 + 1 + 1 + 0 = 26$. La variation totale se décompose en :

$$\text{SST} = \text{SS}_A + \text{SS}_B + \text{SS}_{AB} + \text{SSE}, \qquad \text{SS}_A = 2^k r \, q_A^2$$

Numériquement : $\text{SS}_A = 12 \times 506{,}25 = 6\,075$ (88,7 %), $\text{SS}_B = 675$ (9,9 %), $\text{SS}_{AB} = 75$ (1,1 %) et $\text{SSE} = 26$ (0,4 %), pour un total $\text{SST} = 6\,851$.

On peut aussi donner un intervalle de confiance pour chaque effet. L’écart type de l’erreur vaut $s_e = \sqrt{\text{SSE}/(2^k(r - 1))} = \sqrt{26/8} \approx 1{,}80$ s, et celui d’un effet $s_q = s_e/\sqrt{2^k r} \approx 0{,}52$ s. Avec $t_{0{,}975 \,;\, 8} = 2{,}306$, chaque effet est connu à $\pm 1{,}20$ s près. L’interaction $q_{AB} = 2{,}5 \pm 1{,}2$ s ne contient pas 0 : elle est faible, mais réelle.

Les répétitions sont très groupées, donc l'erreur est petite (±1,2 s sur chaque effet) : même l'interaction de 2,5 s se distingue du bruit.

Piège

Interpréter un petit effet sans répétitions. Sans estimation de l’erreur, rien ne distingue une interaction de 2,5 s d’une fluctuation de mesure. Si la part de l’erreur (SSE/SST) est du même ordre que celle d’un effet, cet effet n’est pas démontré.

Moins d’essais : les plans fractionnaires

Avec 7 facteurs, un plan $2^7$ demande 128 essais. Un plan fractionnaire $2^{k-p}$ n’en réalise qu’une fraction, $2^{k-p}$, choisie astucieusement. Par exemple, un plan $2^{3-1}$ étudie trois facteurs en 4 essais : on réalise le plan $2^2$ en $A$ et $B$, et on fixe le niveau de $C$ par la règle $x_C = x_A x_B$.

Le prix à payer est la confusion (aliasing) : la colonne de $C$ est identique à celle de $AB$, donc l’effet calculé pour $C$ est en réalité $q_C + q_{AB}$. On ne peut pas les séparer. Cela reste acceptable si l’on a de bonnes raisons de penser que les interactions sont faibles devant les effets principaux.

D’où une démarche en deux phases, classique :

1. Criblage : beaucoup de facteurs, deux niveaux chacun, plan $2^k$ ou fractionnaire, pour identifier les quelques facteurs qui comptent ;
2. Approfondissement : ces quelques facteurs seulement, avec plus de niveaux et des répétitions, pour trouver le bon réglage.

En entretien

« Plutôt que de tester les réglages un par un, je fais un plan $2^k$ : chaque facteur à deux niveaux, toutes les combinaisons. La table des signes donne les effets et les interactions, et la part de variation de chacun me dit où investir. Avec des répétitions, j’estime l’erreur pour ne pas interpréter du bruit. »

En résumé

• Faire varier un facteur à la fois coûte peu mais rate les interactions ; le factoriel complet voit tout mais explose.
• Le plan $2^k$ teste $k$ facteurs à deux niveaux codés $\pm 1$ en $2^k$ essais, avec le modèle $y = q_0 + \sum q_A x_A + \sum q_{AB} x_A x_B + \dots$
• Chaque effet se calcule par la table des signes : $q = \frac{1}{2^k} \sum_j s_j y_j$ ; sa part de variation est $2^k q^2/\text{SST}$.
• Le plan $2^k r$ répète chaque essai $r$ fois, ce qui permet d’estimer l’erreur (SSE) et de donner un intervalle de confiance pour chaque effet.
• Les plans fractionnaires $2^{k-p}$ économisent des essais au prix d’effets confondus.
• Démarche en deux phases : criblage de nombreux facteurs, puis approfondissement des facteurs importants.

Exercices

Exercice · Calculer les effets d'un plan 2^2

facile

On mesure le débit (en requêtes/s) d’un serveur selon la taille de la mémoire ($A$ : 8 Go ou 32 Go) et le nombre de fils d’exécution ($B$ : 4 ou 16). Résultats : (8 Go, 4 fils) : 200 ; (32 Go, 4 fils) : 260 ; (8 Go, 16 fils) : 300 ; (32 Go, 16 fils) : 440.

Calculez $q_0$, $q_A$, $q_B$, $q_{AB}$ et la part de variation expliquée par chacun.

Voir la solution

Avec les signes ($A$, $B$) = (−,−), (+,−), (−,+), (+,+) :

• $q_0 = (200 + 260 + 300 + 440)/4 = 300$ ;
• $q_A = (-200 + 260 - 300 + 440)/4 = 50$ ;
• $q_B = (-200 - 260 + 300 + 440)/4 = 70$ ;
• $q_{AB} = (200 - 260 - 300 + 440)/4 = 20$.

$\text{SS}_A = 4 \times 2\,500 = 10\,000$, $\text{SS}_B = 4 \times 4\,900 = 19\,600$, $\text{SS}_{AB} = 4 \times 400 = 1\,600$, donc $\text{SST} = 31\,200$. Les fils expliquent 62,8 % de la variation, la mémoire 32,1 %, l’interaction 5,1 %. L’interaction positive indique que la mémoire supplémentaire est plus utile quand il y a beaucoup de fils (+140 req/s avec 16 fils, contre +60 avec 4).

En code

Python

import numpy as np

xA = np.array([-1, 1, -1, 1])             # mémoire : 8 Go (-1) ou 32 Go (+1)
xB = np.array([-1, -1, 1, 1])             # fils : 4 (-1) ou 16 (+1)
y = np.array([200, 260, 300, 440])        # débit (req/s)

signes = {"I": np.ones(4), "A": xA, "B": xB, "AB": xA * xB}   # table des signes
q = {nom: float(col @ y) / 4 for nom, col in signes.items()}
print(q)                                  # {'I': 300.0, 'A': 50.0, 'B': 70.0, 'AB': 20.0}

SS = {nom: 4 * q[nom] ** 2 for nom in ("A", "B", "AB")}
SST = sum(SS.values())
print(SST, {nom: round(100 * v / SST, 1) for nom, v in SS.items()})
# 31200.0 {'A': 32.1, 'B': 62.8, 'AB': 5.1}

R

xA <- c(-1, 1, -1, 1)                     # mémoire : 8 Go (-1) ou 32 Go (+1)
xB <- c(-1, -1, 1, 1)                     # fils : 4 (-1) ou 16 (+1)
y  <- c(200, 260, 300, 440)               # débit (req/s)

# Avec un codage en ±1, les coefficients de la régression sont exactement les effets
q <- coef(lm(y ~ xA * xB))
q                                         # 300  50  70  20

SS <- 4 * q[-1]^2                         # SS_A, SS_B, SS_AB
round(c(SST = sum(SS), 100 * SS / sum(SS)), 1)   # 31200  32.1  62.8  5.1

Exercice · Le piège du un facteur à la fois

moyen

Reprenez les données du plan $2^3$ du chapitre (format $A$, cache $B$, parallélisme $C$). Un collègue part de la configuration (CSV, sans cache, 1 processus) et teste chaque facteur seul.

1. Quelles trois configurations teste-t-il, et quelles durées observe-t-il ?
2. Il décide d’activer les trois améliorations en additionnant les gains observés. Quelle durée prévoit-il ? Quelle est la durée réelle ?

Voir la solution

1. Il teste (Parquet, sans cache, 1 processus) : 45 s ; (CSV, cache, 1 processus) : 100 s ; (CSV, sans cache, 8 processus) : 70 s. La configuration de départ fait 120 s.
2. Gains observés : 75 s pour le format, 20 s pour le cache, 50 s pour le parallélisme, soit 145 s au total. Il prévoit $120 - 145 = -25$ s, ce qui est absurde ! La durée réelle de (Parquet, cache, 8 processus) est 25 s. Les gains ne s’additionnent pas parce que les facteurs interagissent (surtout $A$ et $C$) : une fois la lecture accélérée, il reste bien moins de travail à paralléliser ou à mettre en cache. Seul un plan factoriel permet de le voir.

En code

Python

import itertools

# Plan 2^3 du chapitre, ordre standard (A varie le plus vite, puis B, puis C)
essais = [(a, b, c) for c, b, a in itertools.product([-1, 1], repeat=3)]
duree = dict(zip(essais, [120, 45, 100, 35, 70, 30, 55, 25]))

depart = duree[(-1, -1, -1)]                          # CSV, sans cache, 1 processus
gains = {"format": depart - duree[(1, -1, -1)],       # un seul facteur change à la fois
         "cache": depart - duree[(-1, 1, -1)],
         "parallélisme": depart - duree[(-1, -1, 1)]}
print(gains)                                          # format 75, cache 20, parallélisme 50
print(depart - sum(gains.values()), duree[(1, 1, 1)]) # -25 (prévision absurde) contre 25 réels

R

# Plan 2^3 du chapitre, ordre standard (A varie le plus vite, puis B, puis C)
plan <- expand.grid(A = c(-1, 1), B = c(-1, 1), C = c(-1, 1))
plan$duree <- c(120, 45, 100, 35, 70, 30, 55, 25)
duree <- function(a, b, c) plan$duree[plan$A == a & plan$B == b & plan$C == c]

depart <- duree(-1, -1, -1)                           # CSV, sans cache, 1 processus
gains <- c(format       = depart - duree(1, -1, -1),  # un seul facteur change à la fois
           cache        = depart - duree(-1, 1, -1),
           parallelisme = depart - duree(-1, -1, 1))
gains                                                 # 75  20  50
c(prevision = depart - sum(gains), reel = duree(1, 1, 1))   # -25  25

Exercice · Un plan fractionnaire

difficile

On veut étudier trois facteurs avec seulement 4 essais, grâce au plan $2^{3-1}$ défini par $x_C = x_A x_B$.

1. Écrivez les quatre essais (niveaux de $A$, $B$, $C$).
2. Avec les durées du plan $2^3$ du chapitre, quelles durées obtient-on pour ces quatre essais ? Calculez l’effet estimé de $C$.
3. Comparez à $q_C$ et $q_{AB}$ du plan complet. Qu’observez-vous ?

Voir la solution

1. Les essais sont (−,−,+), (+,−,−), (−,+,−), (+,+,+) : pour chacun, $x_C$ est le produit $x_A x_B$.
2. Dans le plan complet, ces configurations ont donné 70, 45, 100 et 25 s. L’effet estimé de $C$ vaut $(70 - 45 - 100 + 25)/4 = -12{,}5$ s.
3. Dans le plan complet, $q_C = -15$ et $q_{AB} = 2{,}5$, et l’on retrouve bien $q_C + q_{AB} = -12{,}5$. C’est la confusion : dans le plan fractionnaire, la colonne de $C$ est identique à celle de $AB$, et l’on n’estime que leur somme. Ici l’interaction $AB$ est petite, donc l’estimation reste proche de la vérité ; si elle avait été forte, on aurait pu tirer une conclusion fausse sur $C$.

En code

Python

import itertools
import numpy as np

# Plan complet 2^3 du chapitre (A varie le plus vite)
plan = np.array([(a, b, c) for c, b, a in itertools.product([-1, 1], repeat=3)])
y = np.array([120, 45, 100, 35, 70, 30, 55, 25])
A, B, C = plan.T

# Plan fractionnaire 2^(3-1) : on ne garde que les essais où x_C = x_A * x_B
garde = C == A * B
print(plan[garde].tolist(), y[garde].tolist())
# [[1, -1, -1], [-1, 1, -1], [-1, -1, 1], [1, 1, 1]] [45, 100, 70, 25]
q_C_frac = C[garde] @ y[garde] / 4
q_C, q_AB = C @ y / 8, (A * B) @ y / 8               # effets du plan complet
print(q_C_frac, q_C, q_AB, q_C + q_AB)               # -12.5  -15.0  2.5  -12.5 : confusion

R

# Plan complet 2^3 du chapitre (A varie le plus vite)
plan <- expand.grid(A = c(-1, 1), B = c(-1, 1), C = c(-1, 1))
y <- c(120, 45, 100, 35, 70, 30, 55, 25)

# Plan fractionnaire 2^(3-1) : on ne garde que les essais où x_C = x_A * x_B
garde <- with(plan, C == A * B)
cbind(plan[garde, ], y = y[garde])                   # durées 45, 100, 70, 25
q_C_frac <- sum(plan$C[garde] * y[garde]) / 4
q_C  <- sum(plan$C * y) / 8                          # effets du plan complet
q_AB <- sum(plan$A * plan$B * y) / 8
c(q_C_frac, q_C, q_AB, q_C + q_AB)                   # -12.5  -15  2.5  -12.5 : confusion