TP N° 61

Les plans d’expériences

L’objet de ce TP est de présenter les plans d’expériences et de montrer l’apport d’un outil
d’optimisation globale à leur planification puis à leur résolution.

1 – Présenter la problématique générale des plans d’expériences

2 – Proposer des surfaces de réponse dans le cas linéaire et non linéaire

3 – Présenter des techniques de planification optimale

4 – Appliquer les méthodes des plans d’expériences aux essais en Sûreté de Fonctionnement
1 – Problématique générale des plans d’expériences
Un plan d'expériences est une suite ordonnée d'essais menés afin de connaitre l’influence des
paramètres d'entrée (facteurs) sur les résultats (réponses).
Il permet d’identifier les facteurs les plus influents (analyse de la variance), les interactions entre
ces derniers (corrélations) et d’ajuster puis valider un éventuel modèle mathématique du phénomène
étudié (surface de réponse).
La problématique des plans d’expériences recouvre deux aspects :
- la planification que l’on cherche à élaborer pour obtenir un maximum d'informations à partir
d’un nombre minimum d’essais ou pour rendre robuste une surface de réponse,
- l’exploitation des essais réalisés.
La planification des expériences à réaliser se présente sous la forme d’une matrice d’expérience
constituée de n lignes d’essais, définies chacun par une configuration de valeurs de k facteurs.
X1 X2 … Xk Y
1
2
…
n

Parmi les différents plans d’expériences, les plans factoriels assignent à chaque facteur sa valeur la
plus basse (-1) et sa valeur la plus haute (+1) et limitent ainsi à 2k le nombre d’essais. Ce nombre
peut être encore réduit sans trop pénaliser la précision des résultats (plans factoriels fractionnaires).

2 - Surfaces de réponse
Le modèle mathématique utilisé pour décrire le phénomène étudié peut être de type linéaire ou non
linéaire.

2.1 – Modèle linéaire

Un modèle linéaire est de la forme Y = X β + ε avec Y la matrice des réponses, X la matrice des
facteurs, β la matrice des paramètres à estimer et ε une matrice d’erreur ou de bruit.
Un modèle polynomial est linéarisable en considérant des facteurs multiplicatifs additionnels,
comme le modèle complet suivant du second degré à 3 facteurs :
Y = a0 + a1 X1 + a2 X2 + a3 X3 + a4 X12 + a5 X22 + a6 X32 + a7 X1 X2 + a8 X1 X3 + a9 X2 X3 + ε

La résolution d’un modèle linéaire peut s’effectuer de manière analytique par résolution d'un
système d'équations linéaires : β = (XTX)-1 XTY. Cette estimation est de type moindres carrés :
ε2 = (Xβ-Y)T(Xβ-Y)
dε2/dβ = 2XTXβ - 2XTY
dε2/dβ = 0 ⇒ β = (XTX)-1 XTY

La matrice d’information de Fisher est alors 1/σ2 XTX en considérant un bruit de mesure uniforme.
La matrice de variance-covariance des paramètres (var β) est égale à σ2 (XTX)-1.
Si on considère une erreur de mesure non uniforme, on obtient :

2
β = (XTW-1X)-1 XTW-1Y et var β = (XTWX)-1 avec W la matrice des poids inverse de la matrice de
variance covariance des bruits de mesure.
Traité par l’outil Gencab, l’exemple suivant montre la résolution d’un modèle polynomial du
second degré à 3 facteurs avec estimation des intervalles de confiance.
Plan d'expériences - Modèle linéaire : Y = X P + ε Planification Exploita

min : 0 0 0 Taux de confiance : 60%

max : 1 1 1 Bruit (sigma) : 0,1 Déterminant (Fisher) : 0,0002772
2 2 2
Y X0 X1 X2 X3 x1 x2 x3 x1X2 x1X3 x2X3
2 2 2
1 0,64 1 0,57 0,6 0,2 0,33 0,36 0,04 0,34 0,11 0,12 Y = a0 + a1 X1 + a2 X2 + a3 X3 + a4 X1 + a5 X2 + a6 X3 + a7 X1X2 + a8 X1X3 + a9 X2X3
2 0,82 1 0,46 0,06 0,18 0,21 0 0,03 0,03 0,08 0,01
3 0,59 1 0,64 0,98 0,53 0,41 0,96 0,28 0,63 0,34 0,52 P min max
4 0,94 1 0,47 0,19 0,97 0,22 0,04 0,94 0,09 0,46 0,18 a0 : 0,57 0,53 0,6
5 0,24 1 0,4 0,75 0,78 0,16 0,56 0,61 0,3 0,31 0,58 a1 : 0,67 0,58 0,76
6 0,51 1 0,4 0,15 0,03 0,16 0,02 0 0,06 0,01 0 a2 : -0,8 -0,9 -0,7
7 0,6 1 0,53 0,2 0,99 0,28 0,04 0,98 0,1 0,52 0,19 a3 : -0,8 -0,8 -0,7
8 0,6 1 0,55 0,05 0,84 0,3 0 0,71 0,03 0,46 0,04 a4 : 0,61 0,54 0,67
9 0,81 1 0,73 0,63 0,87 0,54 0,4 0,75 0,46 0,63 0,55 a5 : 0,97 0,92 1,02
10 0,15 1 0,68 0,45 0,77 0,46 0,21 0,6 0,31 0,52 0,35 a6 : 0,74 0,69 0,78
11 0,94 1 0,95 0,96 0,22 0,91 0,93 0,05 0,92 0,21 0,21 a7 : -0,9 -0,9 -0,8
12 0,44 1 0,27 0,47 0,56 0,07 0,22 0,31 0,13 0,15 0,26 a8 : -0,5 -0,6 -0,5
13 0,9 1 0,75 0,12 0,79 0,56 0,01 0,62 0,09 0,59 0,09 a9 : 0,38 0,35 0,42
14 0,15 1 0,5 0,54 0,97 0,25 0,29 0,94 0,27 0,48 0,52
15 0,76 1 0,51 0,8 0,71 0,26 0,64 0,51 0,41 0,36 0,57 x y Var(y) min max
16 0,01 1 0,08 0,45 0,42 0,01 0,21 0,17 0,04 0,03 0,19 1 0,5 0,5 0,5 0,25 0,25 0,25 0,25 0,25 0,25 0,43 0 0,42 0,44
17 0,2 1 0,22 0,26 0,57 0,05 0,07 0,33 0,06 0,13 0,15
18 0,75 1 0,23 0,31 0,2 0,05 0,09 0,04 0,07 0,05 0,06 Var(P)
19 0,65 1 0,45 0,84 0,17 0,2 0,71 0,03 0,38 0,08 0,14 0,05 -0,1 -0,1 -0,1 0,04 0,05 0,02 0,04 0,07 0,02
20 0,66 1 0,25 0,4 0,03 0,06 0,16 0 0,1 0,01 0,01 -0,1 0,33 0,1 0,11 -0,2 -0,1 -0 -0,1 -0,2 -0
21 0,93 1 0,29 0,42 0,18 0,08 0,18 0,03 0,12 0,05 0,08 -0,1 0,1 0,19 0,05 -0 -0,1 -0 -0,1 -0,1 -0
22 0,3 1 0,35 0,97 0,18 0,12 0,95 0,03 0,34 0,06 0,17 -0,1 0,11 0,05 0,13 -0 -0 -0,1 -0 -0,1 -0
23 0,09 1 0,33 0,3 0,33 0,11 0,09 0,11 0,1 0,11 0,1 0,04 -0,2 -0 -0 0,15 0,01 0,01 -0 0,07 0,01
24 0,85 1 0,04 0,88 0,86 0 0,78 0,75 0,04 0,04 0,76 0,05 -0,1 -0,1 -0 0,01 0,1 0,01 0,03 0,07 0,02
25 0,28 1 0,24 0,3 0,1 0,06 0,09 0,01 0,07 0,02 0,03 0,02 -0 -0 -0,1 0,01 0,01 0,08 0,02 -0 -0
26 0,38 1 0,13 0,16 0,66 0,02 0,02 0,44 0,02 0,09 0,1 0,04 -0,1 -0,1 -0 -0 0,03 0,02 0,13 0,04 0
27 0,95 1 0,43 0,87 0,61 0,18 0,76 0,37 0,37 0,26 0,53 0,07 -0,2 -0,1 -0,1 0,07 0,07 -0 0,04 0,21 0,01
28 0,61 1 0,47 0,76 0,96 0,23 0,57 0,93 0,36 0,46 0,73 0,02 -0 -0 -0 0,01 0,02 -0 0 0,01 0,05

Modèle linéaire
Feuille de calcul disponible par double clic sur l’îcone :

2.2 – Modèle non linéaire

La résolution d’un modèle non linéaire impose l’usage d’un outil d’optimisation globale comme
l’exemple de pharmacocinétique suivant, traité par la méthode des moindres carrés.

+ε
- a1 t
Modèle non linéaire Y = a0 e
a0 a1
Σ(Y(ti) - Yi) : 145,6102
2
w: 15 12,51 -0,53

Fisher
Yi ti Y(t) (Y(ti) - Yi)2 I11 I22 I10=I01 D a0 a1
1 25 2 36 120,9977 0,037196 22,456 -0,914 1E-16 Ajustement
2 59 3 61,1 4,306786 0,144804 168,626 -4,88 0,6041
3 180 5 176 17,72513 1,04539 3566,74 -60,2 104,61 1000
4 860 8 858 2,580564 22,85008 214188 -2203 40039 12,29 -0,53 900
800
700
600
500
Y

400
300
200
100
0
0 2 4 6 8 10
20 1,69E+09 1E+14 -3E+10 2E+23 12,29 -0,53 t

3
 Non linéaire 1

L’ajustement peut être également réalisé par la méthode du maximum de vraisemblance, comme
traité ci-après, en considérant que le bruit de mesure est distribué suivant une loi normale. L’outil
Gencab permet alors d’obtenir la matrice de Fisher au moyen d’une méthode numérique qui donne
des résultats voisins de celle calculée de manière analytique.

+ε
- a1 t
Modèle non linéaire Y = a0 e
a0 a1 Taux de confiance : 90%
Σ(Y(ti) - Yi) : 139,6805
2
w: 15 12,29 -0,53 ∑Ln(L) : -3,9862
Min Max
Fisher Vraisemblance a_0 : 12,28531 8,516468 16,05416
2
Yi ti Y(t) (Y(ti) - Yi) I11 I22 I10=I01 D a0 a1 L Ln(L) a_1 : -0,531139 -0,570063 -0,492215
1 25 2 35,5 111,1096 0,037196 22,456 -0,914 1E-16 0,31166 -1,1658
2 59 3 60,5 2,103437 0,144804 168,626 -4,88 0,6041 0,39708 -0,9236 Matrice de Fisher : 22,85008 -2203,224
3 180 5 175 26,21422 1,04539 3566,74 -60,2 104,61 0,37637 -0,9772 -2203,224 214222,3
4 860 8 861 0,253297 22,85008 214188 -2203 40039 12,29 -0,53 0,39872 -0,9195
Variance-covariance : 5,250037 0,053995
0,053995 0,00056

Ajustement
1000
800
600

Y
400
200
20 1,69E+09 1E+14 -3E+10 2E+23 12,29 -0,53 0
0 2 4 6 8 10
t

Non linéaire 2

L’outil permet ainsi d’obtenir un intervalle de confiance pour chacun des paramètres et pour
l’estimation de la réponse d’une configuration de facteurs.

Plan d'expériences - Modèle non linéaire Planificati Exploitati

min : -20 Somme : -3,99 Taux de confiance : 60%

max : 20 Bruit (sigma) : 15
Y X1 y(X) R ln(L) min max
1 25 2 35,54 0,312 -1,17 Va_1 : 12,29 11,32 13,25
2 59 3 60,45 0,397 -0,92 Va_2 : -0,53 -0,54 -0,52
3 180 5 174,9 0,376 -0,98
4 860 8 860,5 0,399 -0,92 x y Var(y) min max
2 35,54 24,57 33,46 37,63 h : 1E-10

35,54 dY/dVa_1 2,893

35,54 dY/dVa_2 -71,1

Matrice de Fisher : 22,85 -2203

-2203 2E+05

Variance-Covariance : 5,252 0,054

0,054 6E-04

4
 Non linéaire 3

3 – Planification optimale
L’optimisation d’un plan d’expériences peut s’effectuer selon divers critères tels que :
- la minimisation du nombre d’essais nécessaires à l’estimation des paramètres du modèle de
surface de réponse,
- la précision de cette estimation,
- La robustesse du modèle à la prédiction dans toutes les configurations de valeurs prises par
les différents facteurs.
Cherchant à minimiser la variance des estimations, la méthode D-optimale répond aux 2 premiers
critères. Des techniques d’occupation optimale de l’espace des configurations d’essais, telle que la
maximisation des distances minimales, répondent au troisième critère.

3.1 – Méthode D-optimale

La méthode D-optimale1 consiste à maximiser le déterminant de la matrice de variance covariance
(inverse de la matrice de Fisher) afin de minimiser la variance des estimateurs et donc les intervalles
de confiance des estimations.
Cette méthode est mise en œuvre dans les exemples précédents afin de déterminer des conditions
optimales d’essais supplémentaires.

3.2 – Occupation optimale de l’espace

La normalisation préalable des différents facteurs permet de définir une notion de distance entre les
conditions d’essais. La maximisation des distances minimales permet alors d’occuper l’espace de
manière homogène comme le montre l’exemple suivant limité à 2 facteurs.

Plan d'expérience - Maximisation de la distance minimale

2
x y min di
10 5 0

2
xi yi di
10 10 25
-10 -10 625 Maximisation de la distance minimale
-10 10 425
10
10 -10 225
0 0 125 8
0 -10 325 6
-10 0 425
4
0 10 125
10 0 25 2
-5 5 225
0
x

5 5 25 -10 -8 -6 -4 -2 0 2 4 6 8 10
5 -5 125 -2
-5 -5 325 -4
10 5 0
-6

-8

-10
y

1
Optimal design, D-optimality, [Link]

5
 Maximisation de la
distance maximale

Cette optimisation des distances, réalisée ici de manière séquentielle, peut l’être également de
manière simultanée en multipliant le nombre de variables (26 dans l’exemple ci-desous).

Plan d'expérience - Maximisation simultanée de la distance minimale

Maximisation de la distance
xi 0 -7,35 4,593 8,376 0,574 -6,31 9,91 -3,57 -10 -1,49 -5,02 -9,71 6,018 5,45
minimale
xi yi 0 -9,79 9,958 1,7 5,774 -3,43 7,62 9,868 1,229 -8,15 4,21 7,991 -9,64 -3,34 10
0 0 149,8 120,3 73,04 33,67 51,63 156,3 110,1 101,5 68,71 42,9 158,1 129,1 40,86 8
-7,35 -9,79 149,8 532,6 379,3 305 41,52 601 400,7 128,4 37,01 201,4 321,6 178,7 205,4
6
4,593 9,958 120,3 532,6 82,5 33,66 298,2 33,73 66,64 289,1 365,1 125,4 208,3 386 177,6
4
8,376 1,7 73,04 379,3 82,5 77,46 242,1 37,4 209,4 337,8 194,4 185,7 366,5 134,1 33,96
0,574 5,774 33,67 305 33,66 77,46 132,1 90,57 33,94 132,4 198,2 33,71 110,6 267,1 106,8 2
-6,31 -3,43 51,63 41,52 298,2 242,1 132,1 385,3 184,3 35,26 45,61 60,02 141,9 190,6 138,4 0

x
9,91 7,62 156,3 601 33,73 37,4 90,57 385,3 186,8 437,2 378,8 234,5 384,9 312,9 140 -10 -5 -2 0 5 10
-3,57 9,868 110,1 400,7 66,64 209,4 33,94 184,3 186,8 116 329,1 34,11 41,18 472,4 255,8
-4
-10 1,229 101,5 128,4 289,1 337,8 132,4 35,26 437,2 116 160,4 33,7 45,81 374,6 259,5
-1,49 -8,15 68,71 37,01 365,1 194,4 198,2 45,61 378,8 329,1 160,4 165,3 328,2 58,57 71,34 -6
-5,02 4,21 42,9 201,4 125,4 185,7 33,71 60,02 234,5 34,11 33,7 165,3 36,28 313,5 166,6 -8
-9,71 7,991 158,1 321,6 208,3 366,5 110,6 141,9 384,9 41,18 45,81 328,2 36,28 558 358,1 -10
6,018 -9,64 129,1 178,7 386 134,1 267,1 190,6 312,9 472,4 374,6 58,57 313,5 558 39,97 y
5,45 -3,34 40,86 205,4 177,6 33,96 106,8 138,4 140 255,8 259,5 71,34 166,6 358,1 39,97
2 2
di max : 33,67 37,01 33,66 33,96 33,66 35,26 33,73 33,94 33,7 37,01 33,7 36,28 39,97 33,96 max (di min) : 33,6625369

Maximisation
simultanée

3 – Essais en Sûreté de fonctionnement

Outre l’application de la méthode D-optimale aux systèmes monocoup traitée dans le TP précédent
(n°60), des essais de fiabilité (exponentielle) et de dégradation (processus de Wiener) accélérés en
température (Arrhenius) sont traités ci-après.

D-optimale - Loi exponentielle accélérée

lbd Ea Température FA E(T) I00 I11 I01
0,44701302 1,36558 -100 0,0853455 16649,079 2019024 791,741 -4,00E+04
Taux de confiance : 90%
k: 0,00008617 Last_lbd Last_Ea
Ln vraisemblance : -300 Déterminant : 7,209E+09 0,0007 0,08746 Min Max
T référence °C Lbd : 0,447013 -0,00064 0,002046
25 Résultats Ea : 1,36558 -0,20031 0,37524
Temps Facteur Densité de
Température Ln(L) FA I00 I11 I01 = I10 lbd Ea Température
cumulé d'accélération probabilité
1 125 627 635201,1336 0 -100 2,3533 2019023,6 95,739056 13903,216 Matrice de Fisher : 6057681 24505,45
2 75 800 2081,28446 0 -100 1,6313 4038047,2 127,04572 21853,618 (numérique) 24505,45 131,8027
3 40 1259 12,78843088 0 -100 1,1773 6057070,8 130,52869 24505,446 0,0007 0,08746 0,75324832
-100 Variance-covariance : 6,66E-07 -0,00012
-0,00012 0,030609

Matrice de Fisher (analytique)

I00 1/λ²
I11 (ln(Fa)/Ea)²
I02 ln(Fa)/λEa

Déterminant
D=I00(I11×I22-I12²)-I11×I02² (I01=I10=0)

Exponentielle

6
 D-optimale - Processus de Wiener linéaire accéléré
m s Ea Temp Fa E(Z) I00 I11 I22 I02 I12 Taux de confiance : 90%
0,1801 0,2944 0,2369 250 52,97 1144 61096 23,068 556289 184344 56,898
Min Max
T référence °C Ln vraisemblance : -9,037 Déterminant : 3E+11 last m last s Last Ea m: 0,1801 0,16518815 0,19501007
25 0,1801 0,2944 0,2369 s: 0,2944 0,09673001 0,49215322
t T°C Z Fa f(t) ln(f(t) Last Fa I00 I11 I22 I02 I12 Ea : 0,2369 0,23694659 0,23694659
1 0 0
2 100 25 20 1 0,1078 -2,227 1 1153,4 23,068 0 0 0 Matrice de Fisher : 12168,1185 0,00207958 13277,2238
3 200 75 80 3,7641 0,0275 -3,592 3,7651 5496,1 46,136 4401,8 4376,1 19,002 (numérique) 0,00207958 69,2038499 44,1604284
4 300 95 190 5,7826 0,04 -3,218 5,7847 12168 69,204 16273 13277 44,16 0,1801 0,2944 0,2369 13277,2238 44,1604284 16487,8973

Variance-covariance : 0,00068585 0,00035301 -0,00055324

0,00035301 0,0146565 -0,00032353
-0,00055324 -0,00032353 0,00050702

Matrice de Fisher (analytique)

I00 nδtα/σ²
I11 2n/σ²
I22 (-1/(2Ea²))×[2Σ(lnαi)-((2δtμ²)/σ²)Σ(αi×(lnαi)²)-Σ(lnαi)²]
I02 [(μδt)/(σ²Ea)]Σ(αi×lnαi)
I12 [1/(Eaσ)]Σlnαi

Déterminant
D=I00(I11×I22-I12²)-I11×I02² (I01=I10=0)

Wiener

La méthode D-optimale à tendance à proposer ici les limites du domaine de variation en

température afin de maximiser la précision des estimations. Il n’en serait pas forcément de même si
l’essai de fiabilité était censuré, ce qui modifierait l’expression de la vraisemblance.
Les matrices calculées de manières analytique et numérique apparaissent toujours voisines.