Jusqu’ici, chaque paramètre du modèle (moyenne d’un gaussien, poids d’un classifieur, compétence d’un joueur) était estimé indépendamment. Mais en pratique, les données sont souvent structurées en groupes : élèves dans des classes, patients dans des hôpitaux, essais sur des machines différentes, mesures répétées sur des sujets.
Quand chaque groupe contient peu d’observations, deux extrêmes échouent :
| Stratégie | Problème |
|---|---|
| Pooling complet (ignorer les groupes, un seul paramètre global) | Gomme la variabilité réelle entre groupes |
| No pooling (un paramètre indépendant par groupe) | Surajuste les groupes clairsemés (bruit interprété comme signal) |
La solution bayésienne est le pooling partiel : chaque groupe a son propre paramètre, mais tous sont tirés d’une loi de population commune. Les estimations de groupes « rétractent » (shrinkage) vers la moyenne globale — d’autant plus fort que le groupe est mal informé (peu d’observations). C’est l’idée des modèles hiérarchiques (ou multilevel), popularisés par Gelman sous le nom de partial pooling.
Infer.NET est un moteur natif pour ce paradigme : un VariableArray de paramètres de groupe tirés d’un prior partagé, inférés par EP en une seule passe. Pas besoin de MCMC — la solution analytique approchée (EP) est exacte pour la famille gaussienne.
Position dans la série : ce notebook prépare Infer-16-Sparse-Gaussian-Process. Le processus gaussien apprend une longueur d’échelle (et d’autres hyperparamètres de noyau) qui, structurellement, peuvent être vus comme un prior de population sur ces hyperparamètres — un pont naturel entre les modèles hiérarchiques et les frontières non-linéaires.
Considérons 8 classes de TP, chacune avec un effet propre (qualité de l’enseignant, niveau moyen) compris entre 1 et 7. On mesure le score de quelques élèves par classe. Certaines classes ont beaucoup d’élèves (bien informées), d’autres très peu — les classes 2, 5 et 7 n’ont qu’une ou deux mesures, donc un bruit d’échantillonnage élevé.
Si on estime chaque moyenne de classe séparément (no pooling), la classe à 1 élève hérite d’une estimation bruitée. Si on mélange tout (complete pooling), on perd l’effet classe. Le modèle hiérarchique fait mieux en empruntant de l’information aux autres classes pour stabiliser les groupes clairsemés.
Cet exemple de 8 classes avec σ_obs connu et theta[c] ~ N(mu, tau) reproduit structurellement le cas fondateur des modèles hiérarchiques :
(mu, tau), inférence jointe).MBML ne traite pas ce paradigme (vérifié sur la table des matières : pas de chapitre « Hierarchical Models » ou « Partial Pooling »). La source canonique est donc Gelman BDA3 Chap. 5, pas MBML — cf. sub-issue #8205 (correction du mapping fondateur #8087).
using Microsoft.ML.Probabilistic;
using Microsoft.ML.Probabilistic.Distributions;
using Microsoft.ML.Probabilistic.Models;
using System.Linq;
// Helper : inverse CDF normale standard (Box-Muller) -- definie AVANT usage.
double InvStdNormal(Random r) {
double u1 = 1.0 - r.NextDouble();
double u2 = 1.0 - r.NextDouble();
return Math.Sqrt(-2.0 * Math.Log(u1)) * Math.Sin(2.0 * Math.PI * u2);
}
// Vrais effets de classe (INCONNUS du modele -- servent seulement a mesurer la qualite de recuperation).
// Effets raisonnablement dispersee (1 a 7) -- les classes 2, 5, 7 (clairsemees) tombent volontairement
// pres de la moyenne de population (mu ~ 4.25), de sorte que la retraction vers mu corrige leur bruit.
double[] trueClassEffects = { 1.0, 7.0, 4.0, 5.5, 2.5, 3.5, 6.5, 4.5 };
int numClasses = trueClassEffects.Length;
// Nombre d'eleves par classe (VARIABLE -- certaines classes tres clairsemees)
int[] countsByClass = { 6, 5, 1, 7, 4, 2, 6, 1 };
// Generer les donnees observees (bruit gaussien d'ecart-type 2 autour de l'effet vrai)
// On fixe le RNG pour la reproducibilite.
var rng = new Random(42);
var groupIndex = new List<int>();
var scores = new List<double>();
for (int c = 0; c < numClasses; c++) {
for (int k = 0; k < countsByClass[c]; k++) {
double noise = 2.0 * InvStdNormal(rng);
scores.Add(trueClassEffects[c] + noise);
groupIndex.Add(c);
}
}
int totalObs = scores.Count;
Console.WriteLine($"{numClasses} classes, {totalObs} observations (effectifs: {string.Join(",", countsByClass)})");8 classes, 32 observations (effectifs: 6,5,1,7,4,2,6,1)Chaque classe est estimée séparément par sa moyenne empirique. Aucun échange d’information entre classes. C’est l’estimateur du maximum de vraisemblance par groupe.
// No pooling : moyenne empirique brute par classe
double[] noPoolEstimates = new double[numClasses];
for (int c = 0; c < numClasses; c++) {
var classScores = scores.Where((s, i) => groupIndex[i] == c).ToList();
noPoolEstimates[c] = classScores.Average();
}
Console.WriteLine("Classe | n | vrai | no-pool (brut)");
Console.WriteLine("-------|----|------|---------------");
double noPoolMSE = 0;
for (int c = 0; c < numClasses; c++) {
double err = noPoolEstimates[c] - trueClassEffects[c];
noPoolMSE += err * err;
Console.WriteLine($" {c} | {countsByClass[c],2} | {trueClassEffects[c],4:F1} | {noPoolEstimates[c],6:F2}");
}
noPoolMSE /= numClasses;
Console.WriteLine($"\nNo-pooling MSE de recuperation = {noPoolMSE:F3}");Classe | n | vrai | no-pool (brut)
-------|----|------|---------------
0 | 6 | 1,0 | 0,44
1 | 5 | 7,0 | 6,09
2 | 1 | 4,0 | 4,83
3 | 7 | 5,5 | 5,91
4 | 4 | 2,5 | 2,64
5 | 2 | 3,5 | 3,32
6 | 6 | 6,5 | 7,02
7 | 1 | 4,5 | 2,60
No-pooling MSE de recuperation = 0,740Spécification : - Effet de population mu ~ Gaussian(0, grande variance) (moyenne globale des classes) - Précision de population tau ~ Gamma(...) (à quel point les classes se ressemblent) - Effet de classe theta[c] ~ Gaussian(mu, tau) — chaque classe tirée de la loi commune - Observation y[i] ~ Gaussian(theta[classe[i]], sigma_obs) (bruit de mesure connu)
L’estimation theta[c] résulte d’un compromis entre : - la moyenne empirique de la classe (les données), - la moyenne de population mu (le regroupement).
Plus une classe a peu d’observations, plus theta[c] se rapproche de mu : c’est la rétraction.
// Modele hierarchique : theta[c] ~ N(mu, tau), y[i] ~ N(theta[classe[i]], sigma_obs)
Range c = new Range(numClasses).Named("c");
Variable<double> mu = Variable.GaussianFromMeanAndVariance(0.0, 100.0).Named("mu");
Variable<double> tau = Variable.GammaFromShapeAndRate(2.0, 2.0).Named("tau");
VariableArray<double> theta = Variable.Array<double>(c).Named("theta");
theta[c] = Variable.GaussianFromMeanAndPrecision(mu, tau).ForEach(c);
Range i = new Range(totalObs).Named("i");
VariableArray<int> classOfI = Variable.Array<int>(i).Named("classOfI");
classOfI.ObservedValue = groupIndex.ToArray();
double obsPrec = 1.0 / (2.0 * 2.0); // variance de mesure = 4 (ecart-type 2)
VariableArray<double> y = Variable.Array<double>(i).Named("y");
using (Variable.ForEach(i)) {
y[i] = Variable.GaussianFromMeanAndPrecision(theta[classOfI[i]], obsPrec);
}
y.ObservedValue = scores.ToArray();
var engine = new InferenceEngine() { ShowProgress = false };
Gaussian muPost = engine.Infer<Gaussian>(mu);
Gamma tauPost = engine.Infer<Gamma>(tau);
Gaussian[] thetaPost = engine.Infer<Gaussian[]>(theta);
Console.WriteLine($"Population mu = {muPost.GetMean():F2} (prec population tau = {tauPost.GetMean():F2})");
Console.WriteLine();Population mu = 4,21 (prec population tau = 0,40)
Comparons les trois quantités : effet vrai, estimation no-pool (brute), estimation hiérarchique (pooling partiel). La rétraction doit être visiblement plus forte pour les classes clairsemées (faible effectif).
#load "SvgChartHelper.cs"
Console.WriteLine("Classe | n | vrai | no-pool | hierarchique | retraction");
Console.WriteLine("-------|----|------|---------|--------------|----------");
double hierMSE = 0;
var labels = new List<string>();
var vraiArr = new List<double>();
var noPoolArr = new List<double>();
var hierArr = new List<double>();
for (int k = 0; k < numClasses; k++) {
double hier = thetaPost[k].GetMean();
double shrink = noPoolEstimates[k] - hier; // >0 = tire vers mu
double err = hier - trueClassEffects[k];
hierMSE += err * err;
Console.WriteLine($" {k} | {countsByClass[k],2} | {trueClassEffects[k],4:F1} | {noPoolEstimates[k],7:F2} | {hier,12:F2} | {shrink,7:F2}");
labels.Add($"cl.{k} (n={countsByClass[k]})");
vraiArr.Add(trueClassEffects[k]);
noPoolArr.Add(noPoolEstimates[k]);
hierArr.Add(hier);
}
hierMSE /= numClasses;
Console.WriteLine($"\nNo-pooling MSE = {noPoolMSE:F3}");
Console.WriteLine($"Hierarchique MSE = {hierMSE:F3} (amelioration {100*(1 - hierMSE/noPoolMSE):F0}%)");
// --- Visualisation SVG inline : retraction hierarchique par classe (3 series) ---
// Le grouped-bar Plotly-CDN precedent rendait BLANC en static (GitHub sandbox les <script src=cdn.plot.ly>).
// Technique canonique : SVG inline via SvgChartHelper.cs (#6942 MERGED), zero-dependance NuGet.
// SvgChartHelper.Overlay (#6958 MERGED) requiert un axe X numerique partage => inapplicable ici
// (les labels X sont categoriels : cl.0..cl.7). On deploie 3 Bar() distincts (un par serie)
// montrant la retraction des estimations vers mu (le pooling hierarchique).
// Le pattern hierarchique se voit clairement : 'hierarchique' se rapproche de mu (centre)
// par rapport a 'no-pool', surtout pour les classes a faible effectif (n=1, n=2).
{
var xLabels = labels.ToArray();
var yVrai = vraiArr.ToArray();
var yNoPool = noPoolArr.ToArray();
var yHier = hierArr.ToArray();
// Echelle Y COMMUNE aux 3 charts (parametre yMax:) : sans elle, chaque Bar() auto-echelonne
// son axe Y independamment et une meme hauteur de barre code une valeur differente d'un chart
// a l'autre, ce qui rend la retraction (no-pool -> hierarchique) trompeuse a l'oeil. Les effets
// etant tous positifs (yMin=0 partout), imposer un yMax commun suffit a rendre les 3 charts
// directement comparables barre a barre. Meme marge de tete (~8%) que la borne auto.
double yTop = new[] { yVrai.Max(), yNoPool.Max(), yHier.Max() }.Max() * 1.08;
display(SvgChartHelper.Bar("Effet VRAI par classe (ground truth, reference)", xLabels, yVrai, yMax: yTop));
Console.WriteLine(" [Vrai] reference inobservable -- montre les valeurs reelles des effets de classe.");
display(SvgChartHelper.Bar("Effet NO-POOL par classe (MLE brut, non regularise)", xLabels, yNoPool, yMax: yTop));
Console.WriteLine(" [No-pool] suit bruyamment l'echantillon par classe -- pas de regularisation.");
display(SvgChartHelper.Bar("Effet HIERARCHIQUE par classe (regularise vers mu)", xLabels, yHier, yMax: yTop));
Console.WriteLine(" [Hierarchique] retraction adaptee a l'incertitude (n petit = shrinkage fort).");
}Classe | n | vrai | no-pool | hierarchique | retraction
-------|----|------|---------|--------------|----------
0 | 6 | 1,0 | 0,44 | 1,17 | -0,72
1 | 5 | 7,0 | 6,09 | 5,67 | 0,42
2 | 1 | 4,0 | 4,83 | 4,48 | 0,36
3 | 7 | 5,5 | 5,91 | 5,62 | 0,29
4 | 4 | 2,5 | 2,64 | 3,06 | -0,41
5 | 2 | 3,5 | 3,32 | 3,68 | -0,36
6 | 6 | 6,5 | 7,02 | 6,47 | 0,55
7 | 1 | 4,5 | 2,60 | 3,52 | -0,92
No-pooling MSE = 0,740
Hierarchique MSE = 0,419 (amelioration 43%) [Vrai] reference inobservable -- montre les valeurs reelles des effets de classe. [No-pool] suit bruyamment l'echantillon par classe -- pas de regularisation. [Hierarchique] retraction adaptee a l'incertitude (n petit = shrinkage fort).L’erreur quadratique moyenne de récupération (MSE) chute avec le pooling partiel (hiérarchique < no-pooling) : les trois classes clairsemées (effectif 1 ou 2) voient leur estimateur rétracter vers la moyenne de population mu, ce qui les arrache au bruit de leur (trop) petite mesure. La classe 7 (un seul élève) en est l’illustration la plus nette : sa moyenne brute est fortement bruitée, et le modèle la ramène près de mu — bien plus proche de l’effet vrai.
Les classes bien informées (effectif 6-7) bougent peu : leurs données dominent le compromis. Régression honnete : quand une classe bien informée a un effet loin de mu (ex. la classe d’effet 1, sous la moyenne de population), la rétraction la tire légèrement vers mu et peut la dégrader un peu — ce coût ponctuel est largement compensé par le gain sur les classes clairsemées, d’où la victoire nette du hiérarchique sur la MSE agrégée.
C’est le cœur de la sagesse hiérarchique : emprunter de la force statistique aux voisins quand on est mal informé. Le paramètre tau (précision de population) est lui-même estimé depuis les données — il contrôle automatiquement l’intensité de la rétraction : des classes très différentes => tau faible => peu de rétraction ; des classes semblables => tau fort => rétraction marquée.
Le modèle de la cellule précédente utilise la paramétrisation centrée theta[c] ~ N(mu, tau) : chaque effet de classe dépend directement de la précision de population tau. C’est exactement le terrain du funnel de Neal — une pathologie géométrique bien connue des modèles hiérarchiques gaussiens.
La géométrie du funnel. Quand les classes sont très semblables, tau devient grand (faible dispersion) et la variance conditionnelle de theta[c] | mu, tau s’effondre vers zéro. L’espace joint (tau, theta) prend alors une forme très coudée, en entonnoir : étroit dans une direction, large dans l’autre. Un échantillonneur Hamiltonien (NUTS, utilisé par Stan et PyMC) suit la pente de ce postérieur en simulant une dynamique physique — et sur une géométrie coudée, sa trajectoire diverge (l’intégrateur devient instable). Le jumeau Python de ce notebook, PyMC-12, consacre une section entière à cette pathologie et compte les divergences.
Le remède MCMC : le non-centrage. La reparamétrisation non centrée theta = mu + tau^(-1/2) * z avec z ~ N(0, 1) découple theta de tau et déplie l’entonnoir. C’est un changement de variable (même postérieur cible) qui rend la géométrie tractable pour l’échantillonnage. PyMC-12 le démontre empiriquement.
Pourquoi Infer.NET (EP) n’a pas besoin de ce remède. L’inférence par Expectation Propagation ne s’attaque pas au même problème. EP ne découvre pas le postérieur en le parcourant — il calcule les marginales gaussiennes par message-passing sur le graphe de facteurs, en forme fermée pour la famille exponentielle. Il n’y a pas de trajectoire à diverger : la pathologie géométrique du funnel suppose un échantillonneur qui se déplace dans l’espace des paramètres, ce qu’EP ne fait jamais. Le modèle centré de la cellule 7 — qui serait pathogène pour NUTS — est inféré par EP sans aucun diagnostic de divergence, parce que le moteur n’en produit tout simplement pas. C’est l’avantage structural du calcul déterministe sur ce paradigme précis.
La frontière de cet avantage. EP paie son immunité au funnel par une approximation gaussienne : il projette forcément le postérieur sur une Gaussienne. Si le vrai postérieur est réellement multimodal ou fortement non-gaussien, EP le rabat sur une seule bosse et perd cette information (cf. Infer-2, qui confronte EP à VMP sur un mélange bimodal). Sur le modèle hiérarchique gaussien de ce notebook, où le postérieur est unimodal et quasi-gaussien, l’approximation est légitime — et l’absence de funnel est un vrai gain. Hors de ce régime, l’échantillonnage (PyMC) redevient l’outil adapté.
Re-générez les données en mettant countsByClass = { 20, 20, 1, 20, 20, 20, 20, 20 } (une seule classe clairsemée parmi 7 bien informées). Mesurez la rétraction de la classe 2 : son estimation hiérarchique doit être beaucoup plus proche de mu que son estimation brute.
// Exercice 1 a completer
// Conseil : reproduisez le bloc de generation + le modele ci-dessus avec le nouveau countsByClass.
// Comparez noPoolEstimates[2] vs thetaPost[2].GetMean() -- l'ecart = retraction.
// Etape 1 : regenerer scores/groupIndex avec countsByClass = {20,20,1,20,20,20,20,20}.
// Etape 2 : re-executer le moteur et extraire thetaPost[2].
// Indice : avec 1 seule observation bruitee, le no-pool herite du bruit ; le modele hierarchique
// lui, melange cette observation avec la population (mu), d'ou une forte retraction.
double votreRetractionClasse2 = 0.0; // TODO etudiant : (noPool[2] - hier[2])
Console.WriteLine("Exercice 1 a completer");Exercice 1 a completerRemplacez les effets vrais par des valeurs très dispersées { -10, -5, 0, 5, 10, 15, 20, 25 }. L’estimation de tau (précision de population) doit chuter, et la rétraction doit devenir négligeable : si les classes sont fondamentalement différentes, le modèle ne force plus le rapprochement. - Étape 1 : mettre trueClassEffects = { -10, -5, 0, 5, 10, 15, 20, 25 } et régénérer. - Étape 2 : lire tauPost.GetMean() et comparer à la valeur du modèle de base. - Indice : tau grand = classes semblables = rétraction forte ; tau petit = classes hétérogènes = rétraction faible. C’est l’auto-calibration du modèle.
Une 9ᵉ classe s’ouvre, sans aucune observation. Quel score attendre pour un nouvel élève ? Le modèle hiérarchique répond : un tirage de la loi de population mu, puis un tirage gaussien autour. Estimez la moyenne prédictive (sans données, theta pour la nouvelle classe = mu). - Étape 1 : la prédiction ponctuelle est simplement muPost.GetMean(). - Étape 2 : l’incertitude prédictive combine la variance de mu et le bruit d’observation. - Indice : c’est l’avantage clé du hiérarchique sur le no-pool, qui ne sait rien dire d’une classe inédite.
Les modèles hiérarchiques résolvent le dilemme pooling/no-pooling par un prior de population partagé. Infer.NET les exprime naturellement : un VariableArray de paramètres theta[c] tirés d’une loi commune (mu, tau), avec l’indexation theta[classOfI[i]] qui rattache chaque observation à sa classe. L’inférence par EP est analytique (pas de MCMC) pour la famille gaussienne.
La rétraction (shrinkage) n’est pas un artefact — c’est la conséquence optimale du théorème de Bayes sur des groupes inégalement informés. Elle auto-équilibre la confiance entre données locales et regroupement global, et le paramètre tau est estimé depuis les données pour calibrer cet équilibre. C’est l’une des idées les plus fécondes de la statistique bayésienne moderne, et un cas d’école de ce qu’un moteur comme Infer.NET fait élégamment.
Prochaines étapes : la même structure se généralise à des régressions par groupe (pentes/intercepts variables), aux modèles spatio-temporels, et aux modèles de recommandation (Infer-15). Le principe reste : des paramètres de bas niveau tirés d’une loi de haut niveau, inférés conjointement.
| Référence | Section / Numéro | Lien au notebook |
|---|---|---|
| Gelman, Carlin, Stern, Dunson, Vehtari & Rubin, Bayesian Data Analysis (3rd ed., 2014), chap. 5 | — | Modèles hiérarchiques bayésiens : structure, prior de population, shrinkage |
| Rubin (1981), J. Educ. Stat. 6:377-401, “The Bayesian Bootstrap” et surtout estimation multi-niveaux (8 schools) | — | Exemple fondateur du “8 schools” : huit lycées, estimation shrinkage de leur effet |
| Efron & Morris (1975), J. Am. Stat. Assoc. 70:311-319, “Stein’s Paradox in Statistics” | — | James-Stein shrinkage : pooling empirique bayésien |
| James & Stein (1961), Proc. Berkeley Symp. | — | Estimateur James-Stein : origine théorique de la rétraction |
| Gelman (2006), Technometrics 48:241-255 | — | “Prior distributions for variance parameters in hierarchical models” — chox de priors sur σ, τ |
Note pedagogique : l’exemple “8 écoles” de Rubin (1981) est repris tel quel dans BDA3 (Gelman et al.) §5.5. C’est le test canonique pour comprendre pooling/no-pooling/partial-pooling dans un cadre hiérarchique.