polyband

The problemLe problème

InstallInstallation

One command, straight from GitHub. Python 3.9 or newer, with numpy and scipy; matplotlib is only needed for the plotting helpers.

Une seule commande, directement depuis GitHub. Python 3.9 ou plus récent, avec numpy et scipy ; matplotlib n'est requis que pour les fonctions de tracé.

pip install git+https://github.com/eartigau/polyband.git

That installs whatever is on main at the moment you run it. For anything you need to reproduce later, pin a released version instead:

Cette commande installe l'état de main au moment où vous la lancez. Pour tout ce que vous devrez reproduire plus tard, épinglez plutôt une version publiée :

pip install "git+https://github.com/eartigau/polyband.git@v0.2.0"

Quick startDémarrage rapide

order_mean is the degree of the trend, order_width the degree of ln(sigma). They are independent, which is the whole point: a curved trend with a constant-width band is an ordinary combination, and so is a straight trend whose scatter grows.

order_mean est le degré de la tendance, order_width celui de ln(sigma). Ils sont indépendants, et c'est là tout l'intérêt : une tendance courbe avec une enveloppe de largeur constante est une combinaison parfaitement ordinaire, tout comme une tendance droite dont la dispersion augmente.

Copy the whole block and run it: it builds its own data, so it works before you have plugged in yours. Every number printed below is what you should actually see.

Copiez le bloc entier et exécutez-le : il fabrique ses propres données, donc il fonctionne avant même que vous n'y branchiez les vôtres. Chaque valeur indiquée ci-dessous est bien celle que vous devriez obtenir.

import numpy as np
import matplotlib.pyplot as plt
from polyband import fit_polyband, plot_polyband

# ---------------------------------------------------------------------
# 1. The data. Replace this whole block with your own x and y.
# ---------------------------------------------------------------------
rng = np.random.default_rng(0)

# 400 points spread over the x range you care about.
x = rng.uniform(0, 10, 400)

# The curve the points scatter around: a parabola that bends back down.
trend = 3.0 + 1.1 * x - 0.075 * x**2

# How wide that scatter is. The point of the package is that this is NOT
# a constant: it grows from 0.30 at x = 0 to 4.95 at x = 10. A single
# global RMS would be wrong at both ends at once.
sigma = np.exp(-1.2 + 0.28 * x)

# One draw per point, each from its own local sigma.
y = trend + rng.normal(0.0, sigma)

# ---------------------------------------------------------------------
# 2. The fit. Two degrees, chosen independently of each other.
# ---------------------------------------------------------------------
#   order_mean=2  : the trend is a parabola, so degree 2.
#   order_width=1 : sigma is an exponential of x, so ln(sigma) is a
#                   straight line, so degree 1.
# Both polynomials are fitted at the same time by maximum likelihood,
# which is what lets the trend be weighted by the local scatter.
fit = fit_polyband(x, y, order_mean=2, order_width=1)

# ---------------------------------------------------------------------
# 3. Reading the result.
# ---------------------------------------------------------------------
print(fit.summary())          # degrees, N, x range, log L, AIC and BIC

print(fit.predict(5.0))       # the trend at x = 5              -> 6.610
print(fit.scatter(5.0))       # half-width of the 1-sigma band  -> 1.234
print(fit.envelope(5.0))      # the two band edges at x = 5     -> (5.376, 7.844)
print(fit.zscore(5.0, 12.3))  # is y = 12.3 unusual at x = 5?   -> 4.61 sigma

# ---------------------------------------------------------------------
# 4. Check it before you trust it.
# ---------------------------------------------------------------------
# About 68% of the points belong inside the 1-sigma band and 95% inside
# 2 sigma. A band nobody has checked is decoration, not a measurement.
for nsigma, observed, expected in fit.coverage(x, y):
    print(f"{nsigma:.0f} sigma: {observed:.1%} observed vs {expected:.1%} expected")
# 1 sigma: 66.0% observed vs 68.3% expected
# 2 sigma: 96.2% observed vs 95.4% expected
# 3 sigma: 99.8% observed vs 99.7% expected

# ---------------------------------------------------------------------
# 5. Draw it.
# ---------------------------------------------------------------------
# nsigma=(1, 2) nests two bands. plot_polyband draws on the axis you give
# it, never calls show() or legend() itself, and hands back every artist
# it created so you can restyle afterwards.
art = plot_polyband(fit, x, y, nsigma=(1, 2))
art.ax.legend(handles=art.legend_handles)
plt.show()

import numpy as np
import matplotlib.pyplot as plt
from polyband import fit_polyband, plot_polyband

# ---------------------------------------------------------------------
# 1. Les données. Remplacez tout ce bloc par vos propres x et y.
# ---------------------------------------------------------------------
rng = np.random.default_rng(0)

# 400 points répartis sur l'intervalle en x qui vous intéresse.
x = rng.uniform(0, 10, 400)

# La courbe autour de laquelle les points se dispersent : une parabole
# qui finit par redescendre.
trend = 3.0 + 1.1 * x - 0.075 * x**2

# L'ampleur de cette dispersion. Tout l'intérêt du module est qu'elle
# n'est PAS constante : elle passe de 0,30 en x = 0 à 4,95 en x = 10.
# Un unique écart-type global se tromperait aux deux bouts à la fois.
sigma = np.exp(-1.2 + 0.28 * x)

# Un tirage par point, chacun selon son sigma local.
y = trend + rng.normal(0.0, sigma)

# ---------------------------------------------------------------------
# 2. L'ajustement. Deux degrés, choisis indépendamment l'un de l'autre.
# ---------------------------------------------------------------------
#   order_mean=2  : la tendance est une parabole, donc degré 2.
#   order_width=1 : sigma est une exponentielle de x, donc ln(sigma) est
#                   une droite, donc degré 1.
# Les deux polynômes sont ajustés simultanément par maximum de
# vraisemblance, ce qui permet de pondérer la tendance par la
# dispersion locale.
fit = fit_polyband(x, y, order_mean=2, order_width=1)

# ---------------------------------------------------------------------
# 3. Lire le résultat.
# ---------------------------------------------------------------------
print(fit.summary())          # degrés, N, intervalle en x, log L, AIC et BIC

print(fit.predict(5.0))       # la tendance en x = 5              -> 6,610
print(fit.scatter(5.0))       # demi-largeur de l'enveloppe à 1s  -> 1,234
print(fit.envelope(5.0))      # les deux bords en x = 5           -> (5,376 ; 7,844)
print(fit.zscore(5.0, 12.3))  # y = 12,3 est-il atypique en x = 5 ? -> 4,61 sigma

# ---------------------------------------------------------------------
# 4. Vérifiez avant de faire confiance.
# ---------------------------------------------------------------------
# Environ 68 % des points doivent tomber dans l'enveloppe à 1 sigma, et
# 95 % dans celle à 2 sigma. Une enveloppe que personne n'a vérifiée est
# une décoration, pas une mesure.
for nsigma, observed, expected in fit.coverage(x, y):
    print(f"{nsigma:.0f} sigma: {observed:.1%} observé vs {expected:.1%} attendu")
# 1 sigma : 66,0 % observé vs 68,3 % attendu
# 2 sigma : 96,2 % observé vs 95,4 % attendu
# 3 sigma : 99,8 % observé vs 99,7 % attendu

# ---------------------------------------------------------------------
# 5. Tracer.
# ---------------------------------------------------------------------
# nsigma=(1, 2) imbrique deux enveloppes. plot_polyband dessine sur l'axe
# que vous lui donnez, n'appelle jamais show() ni legend() de lui-même, et
# renvoie tous les objets graphiques créés pour que vous puissiez les
# restyler ensuite.
art = plot_polyband(fit, x, y, nsigma=(1, 2))
art.ax.legend(handles=art.legend_handles)
plt.show()

The key distinctionLa distinction essentielle

The band is not the uncertainty on the trendL'enveloppe n'est pas l'incertitude sur la tendance

The modelLe modèle

Under the hoodSous le capot

The objective, and the two terms that fight over itLa fonction objectif, et les deux termes qui s'y opposent

Why an underestimate here is not paid for by an overestimate therePourquoi une sous-estimation ici n'est pas compensée par une surestimation ailleurs

Where compensation is genuinely possibleOù la compensation est réellement possible

Finding the minimum in practiceTrouver le minimum en pratique

• Rescaling. x is mapped to t = (x − x₀)/Δ on [−1, 1] before the Vandermonde matrices are built. Raw values in the thousands are hopelessly ill-conditioned by degree 3.
• Starting guess. Ordinary least squares for the trend, then a polynomial fit to ln|residual| for the width. For Gaussian z, E[ln|z|] = ln σ − (γ + ln 2)/2 ≈ ln σ − 0.6352, so that offset is removed; otherwise the optimiser starts from a band a factor of two too narrow.
• Two optimisers. Nelder-Mead first, because the likelihood is mildly non-quadratic in the width coefficients and it does not care; BFGS then polishes from a good point.
• Convergence is judged on the objective, not on the exit flag. BFGS routinely reports precision loss when handed an already-converged start, which is exactly the situation here.
• Covariance comes from an explicit central-difference Hessian, inverted. The inverse Hessian BFGS accumulates is not trustworthy here, since it exits having built almost no curvature information. n_bootstrap resamples instead, at the cost of speed and with no assumption of local quadraticity.
• Bias correction, last. b₀ is raised by ½ ln(N/(N − p)), which undoes the downward bias of maximum-likelihood widths. It is the same correction as dividing by N − 1 rather than N in a plain standard deviation.

• Remise à l'échelle. x est ramené sur t = (x − x₀)/Δ dans [−1, 1] avant la construction des matrices de Vandermonde. Des valeurs brutes de l'ordre du millier sont irrémédiablement mal conditionnées dès le degré 3.
• Point de départ. Moindres carrés ordinaires pour la tendance, puis ajustement polynomial de ln|résidu| pour la largeur. Pour un z gaussien, E[ln|z|] = ln σ − (γ + ln 2)/2 ≈ ln σ − 0,6352 ; ce décalage est donc retiré, sans quoi l'optimiseur partirait d'une enveloppe deux fois trop étroite.
• Deux optimiseurs. Nelder-Mead d'abord, car la vraisemblance est légèrement non quadratique en les coefficients de largeur et cela lui est indifférent ; BFGS peaufine ensuite à partir d'un bon point.
• La convergence est jugée sur l'objectif, non sur le code de sortie. BFGS signale couramment une perte de précision quand on lui donne un point déjà convergé, ce qui est exactement le cas ici.
• La covariance provient d'un hessien numérique par différences centrées, inversé. Le hessien inverse accumulé par BFGS n'est pas fiable ici, puisqu'il s'arrête sans avoir construit d'information de courbure. n_bootstrap rééchantillonne à la place, plus lentement mais sans supposer la vraisemblance localement quadratique.
• Correction de biais, en dernier. b₀ est relevé de ½ ln(N/(N − p)), ce qui annule le biais vers le bas des largeurs par maximum de vraisemblance. C'est la même correction que diviser par N − 1 plutôt que par N dans un écart-type ordinaire.

What the two degrees controlCe que contrôlent les deux degrés

Each row sweeps one degree while the other is held at the value the data were generated with, so each is a clean one-dimensional experiment. Both rows run past the correct answer on purpose: too few degrees is not the only way to be wrong.

Chaque rangée balaie un degré pendant que l'autre est maintenu à la valeur ayant servi à générer les données : chacune est donc une expérience proprement unidimensionnelle. Les deux rangées dépassent volontairement la bonne réponse, car manquer de degrés n'est pas la seule façon de se tromper.

Choosing the degreesChoisir les degrés

If you have no reason to prefer particular degrees, let BIC^[2] scan them. It penalises extra parameters more firmly than AIC^[1], which is what you want here: the failure mode to avoid is a width polynomial flexible enough to chase noise.

Si rien ne justifie un choix particulier de degrés, laissez le BIC^[2] les balayer. Il pénalise les paramètres supplémentaires plus fermement que l'AIC^[1], ce qui est souhaitable ici : le mode de défaillance à éviter est un polynôme de largeur assez souple pour épouser le bruit.

from polyband import select_orders

# Fits every combination from (0, 0) up to the two maxima and scores each one.
# Combinations that fail to fit are skipped rather than raising, so this is
# safe to point at data you have not inspected yet.
fit, table = select_orders(x, y, max_order_mean=4, max_order_width=2)

# `fit` is the winner, already fitted and ready to use. No need to refit.
print(fit.order_mean, fit.order_width)

# `table` is every combination that worked, best first, as
# (order_mean, order_width, criterion). Look at the runners-up rather than
# just the winner: if the top few are within about 2 of each other they are
# a tie, and the simplest of them is the one to keep.
for order_mean, order_width, bic in table[:5]:
    print(f"order_mean={order_mean} order_width={order_width}  BIC={bic:.1f}")

# Pass criterion="aic" to be more permissive, or any fit_polyband argument
# (log_y, nu, yerr) to scan under the same assumptions as your final fit.
fit, table = select_orders(x, y, max_order_mean=4, max_order_width=2, nu=4)

from polyband import select_orders

# Ajuste toutes les combinaisons de (0, 0) jusqu'aux deux maxima et attribue
# un score à chacune. Les combinaisons qui échouent sont ignorées plutôt que
# de lever une exception : on peut donc lancer ceci sur des données non
# encore inspectées.
fit, table = select_orders(x, y, max_order_mean=4, max_order_width=2)

# `fit` est le gagnant, déjà ajusté et prêt à l'emploi. Inutile de réajuster.
print(fit.order_mean, fit.order_width)

# `table` contient toutes les combinaisons ayant fonctionné, meilleure en
# tête, sous la forme (order_mean, order_width, critère). Regardez les
# suivants et pas seulement le gagnant : si les premiers sont à moins de 2
# les uns des autres, c'est une égalité, et c'est le plus simple qu'il faut
# retenir.
for order_mean, order_width, bic in table[:5]:
    print(f"order_mean={order_mean} order_width={order_width}  BIC={bic:.1f}")

# Passez criterion="aic" pour être plus permissif, ou n'importe quel argument
# de fit_polyband (log_y, nu, yerr) pour balayer sous les mêmes hypothèses
# que votre ajustement final.
fit, table = select_orders(x, y, max_order_mean=4, max_order_width=2, nu=4)

Rule of thumb: treat a BIC difference below about 2 as a tie^[3] and keep the simpler model. An over-flexible width is not a harmless extra degree of freedom; it produces a band that wiggles to chase noise and then misstates how unusual any given point is.

Règle empirique : considérez un écart de BIC inférieur à environ 2 comme une égalité^[3] et conservez le modèle le plus simple. Une largeur trop souple n'est pas un degré de liberté anodin : elle produit une enveloppe qui ondule pour suivre le bruit, et donne ensuite une mauvaise mesure du caractère atypique de chaque point.

Quantities that span orders of magnitudeQuantités s'étendant sur plusieurs ordres de grandeur

When the scatter is multiplicative rather than additive, y itself is the wrong variable to fit. Pass log_y=True and polyband works on log10(y) internally while handing results back in the units of y, so nothing downstream has to know about the transform.

Quand la dispersion est multiplicative plutôt qu'additive, y n'est pas la bonne variable à ajuster. Passez log_y=True et polyband travaille en interne sur log10(y) tout en renvoyant les résultats dans les unités de y, de sorte que rien en aval n'a besoin de connaître la transformation.

The one place the transform stays visible is the width: fit.scatter() returns dex, because a band symmetric in log space is asymmetric in linear space and there is no single number to report otherwise. Use fit.envelope() for the edges, which does come back in the units of y.

Le seul endroit où la transformation reste visible est la largeur : fit.scatter() renvoie des dex, car une enveloppe symétrique en espace logarithmique est asymétrique en espace linéaire et aucun nombre unique ne peut la décrire autrement. Utilisez fit.envelope() pour les bords, qui sont bien renvoyés dans les unités de y.

Samples with outliersÉchantillons avec valeurs aberrantes

Under a Gaussian likelihood, a handful of far-flung points drags the band open across the whole x range. Setting nu^[5] switches to a Student-t likelihood, which stops letting those few points set the width. Values of 3 to 5 are strongly resistant; larger values approach the Gaussian case.

Sous une vraisemblance gaussienne, une poignée de points très éloignés élargit l'enveloppe sur tout l'intervalle en x. Définir nu^[5] bascule vers une vraisemblance de Student, qui empêche ces quelques points de dicter la largeur. Des valeurs de 3 à 5 sont très résistantes ; des valeurs plus grandes tendent vers le cas gaussien.

RobustnessRobustesse

One weight, applied to both polynomials at onceUn seul poids, appliqué simultanément aux deux polynômes

And it is already a mixture model. A Student-t is a continuous scale mixture of Gaussians: if y | λ ~ N(μ, σ²/λ) with λ ~ Gamma(ν/2, ν/2), then y ~ t_ν. The posterior of that per-point precision is E[λ_i | y_i] = (ν+1)/(ν+z_i²), which is u_i exactly. So the weights above are the EM weights of a mixture, with the mixing integrated out analytically. "Student-t or a mixture model" is not really a choice between two things.

Et c'est déjà un modèle de mélange. Une loi de Student est un mélange continu de gaussiennes en échelle : si y | λ ~ N(μ, σ²/λ) avec λ ~ Gamma(ν/2, ν/2), alors y ~ t_ν. La loi a posteriori de cette précision par point vaut E[λ_i | y_i] = (ν+1)/(ν+z_i²), soit exactement u_i. Les poids ci-dessus sont donc les poids EM d'un mélange, dont la loi de mélange a été intégrée analytiquement. « Student ou modèle de mélange » n'est pas vraiment un choix entre deux choses.

How deviant they are, and how many there areÀ quel point ils dévient, et combien ils sont

The same two axes, measuredLes deux mêmes axes, mesurés

What nu costs, and the calibration trapCe que coûte nu, et le piège de calibration

Median over 60 clean samples of 600 points. Read the middle column as the deflation to undo if you want a standard deviation, and the right-hand column as the nsigma to ask for if you want a 68.3% interval. Better still, measure it on your own data, which takes one line:

Médiane sur 60 échantillons propres de 600 points. Lisez la colonne du milieu comme la déflation à annuler si vous voulez un écart-type, et celle de droite comme le nsigma à demander si vous voulez un intervalle à 68,3 %. Mieux encore, mesurez-la sur vos propres données, ce qui tient en une ligne :

# With nu set, one fitted width is no longer a 68% interval. Ask coverage()
# what multiple actually is, on your own data rather than from the table.
robust = fit_polyband(x, y, 2, 1, nu=4)

for nsigma, observed, expected in robust.coverage(x, y, nsigma=(1.0, 1.19, 2.0)):
    print(f"{nsigma:.2f} x width: {observed:.1%} inside")

# Then use that multiple everywhere you want a 1-sigma-equivalent band.
lo, hi = robust.envelope(robust.grid(), nsigma=1.19)

# Avec nu défini, une largeur ajustée n'est plus un intervalle à 68 %.
# Demandez à coverage() quel multiple l'est réellement, sur vos propres
# données plutôt que d'après le tableau.
robust = fit_polyband(x, y, 2, 1, nu=4)

for nsigma, observed, expected in robust.coverage(x, y, nsigma=(1.0, 1.19, 2.0)):
    print(f"{nsigma:.2f} x largeur : {observed:.1%} à l'intérieur")

# Utilisez ensuite ce multiple partout où vous voulez une enveloppe
# équivalente à 1 sigma.
lo, hi = robust.envelope(robust.grid(), nsigma=1.19)

Against sigma clipping and against a fitted mixtureFace à l'écrêtage sigma et à un mélange ajusté

Points with measurement errorsPoints avec barres d'erreur

If your points carry uncertainties, the observed spread is the intrinsic spread and the measurement noise added in quadrature. Usually the intrinsic part is the one you want to describe. Pass yerr and the fit separates them inside the likelihood.

Si vos points portent des incertitudes, la dispersion observée est la somme quadratique de la dispersion intrinsèque et du bruit de mesure. C'est généralement la part intrinsèque que vous cherchez à décrire. Passez yerr et l'ajustement les sépare au sein de la vraisemblance.

Compared to binningComparaison avec le binning

Matplotlib integrationIntégration matplotlib

plot_polyband draws on the axis you give it, never touches figure-level state, never calls show() or legend() by itself, and returns every artist it created. That makes it safe to call twice on the same axis, which is how you overlay two fits.

plot_polyband dessine sur l'axe que vous lui donnez, ne touche jamais à l'état global de la figure, n'appelle jamais show() ni legend() de lui-même, et renvoie tous les objets graphiques qu'il a créés. On peut donc l'appeler deux fois sur le même axe, ce qui est la façon de superposer deux ajustements.

fig, ax = plt.subplots(figsize=(9, 5))

art = plot_polyband(
    fit, x, y, ax=ax,                     # ax is yours; nothing global is touched
    nsigma=(1, 2, 3),                     # three nested bands, widest most transparent
    color="#c77dff",                      # trend and bands
    point_kw=dict(s=14, marker="^"),      # forwarded verbatim to ax.scatter
    trend_kw=dict(linewidth=3.0),         # forwarded verbatim to ax.plot
    band_kw=dict(hatch="//"),             # forwarded verbatim to ax.fill_between
    show_trend_error=True,                # add the (much narrower) error on the curve
)

# Everything it drew comes back, so restyle after the fact instead of
# passing more arguments in.
art.trend.set_zorder(10)
art.points.set_alpha(0.15)
ax.legend(handles=art.legend_handles)     # it never calls legend() itself

# Call it twice on the same axis to overlay a second fit. show_points=False
# stops the scatter being drawn again on top of itself.
plot_polyband(robust_fit, ax=ax, show_points=False,
              color="tab:green", label_prefix="Student-t: ")

# Curves stop where the data stop. Ask for an extension deliberately, and
# style it so the reader can see it is an extrapolation.
plot_polyband(fit, ax=ax, show_points=False, extrapolate=0.15,
              trend_kw=dict(linestyle=":"), labels=False)

fig, ax = plt.subplots(figsize=(9, 5))

art = plot_polyband(
    fit, x, y, ax=ax,                     # l'axe est le vôtre, rien de global n'est touché
    nsigma=(1, 2, 3),                     # trois enveloppes imbriquées, la plus large la plus transparente
    color="#c77dff",                      # tendance et enveloppes
    point_kw=dict(s=14, marker="^"),      # transmis tel quel à ax.scatter
    trend_kw=dict(linewidth=3.0),         # transmis tel quel à ax.plot
    band_kw=dict(hatch="//"),             # transmis tel quel à ax.fill_between
    show_trend_error=True,                # ajoute l'erreur (bien plus étroite) sur la courbe
)

# Tout ce qui a été dessiné est renvoyé : restylez après coup plutôt que de
# passer davantage d'arguments.
art.trend.set_zorder(10)
art.points.set_alpha(0.15)
ax.legend(handles=art.legend_handles)     # la fonction n'appelle jamais legend() elle-même

# Appelez-la deux fois sur le même axe pour superposer un second ajustement.
# show_points=False évite de redessiner le nuage par-dessus lui-même.
plot_polyband(robust_fit, ax=ax, show_points=False,
              color="tab:green", label_prefix="Student-t: ")

# Les courbes s'arrêtent où s'arrêtent les données. Demandez une extension
# délibérément, et stylez-la pour que le lecteur voie qu'il s'agit d'une
# extrapolation.
plot_polyband(fit, ax=ax, show_points=False, extrapolate=0.15,
              trend_kw=dict(linestyle=":"), labels=False)

Curves stop exactly where the data stop. A polynomial band has no business being drawn where nothing constrains it, so extrapolate defaults to 0 and any extension is something you have to ask for deliberately.

Les courbes s'arrêtent exactement là où s'arrêtent les données. Une enveloppe polynomiale n'a rien à faire là où rien ne la contraint : extrapolate vaut 0 par défaut, et toute extension doit être demandée délibérément.

Missing and non-finite dataDonnées manquantes et non finies

ReferenceRéférence

Fitting optionsOptions d'ajustement

On the fit objectSur l'objet ajusté

NotesNotes

What AIC, BIC and t₃ actually areCe que sont réellement l'AIC, le BIC et t₃

1. [1]
   AIC, the Akaike information criterion, is 2p − 2 ln L, where p is the number of free parameters and L the likelihood at the optimum. Adding a parameter always raises the likelihood, so the fit alone can never tell you to stop; the 2p term is the price. What Akaike showed is that this particular price makes the criterion an estimate of how much information is lost by using the model instead of the truth, up to a constant common to all the candidates. Two consequences follow. Only differences between AIC values mean anything, never one value on its own. And AIC aims at predicting well, not at identifying a true model: it does not assume the truth is among the candidates, and it does not become certain as the sample grows. See [1a], [1c].
   L'AIC, critère d'information d'Akaike, vaut 2p − 2 ln L, où p est le nombre de paramètres libres et L la vraisemblance à l'optimum. Ajouter un paramètre augmente toujours la vraisemblance : l'ajustement seul ne peut donc jamais vous dire de vous arrêter, et le terme 2p est le prix à payer. Ce qu'Akaike a montré, c'est que ce prix-là fait du critère une estimation de la quantité d'information perdue en utilisant le modèle plutôt que la réalité, à une constante près commune à tous les candidats. Deux conséquences. Seules les différences d'AIC ont un sens, jamais une valeur isolée. Et l'AIC vise à bien prédire, non à identifier un modèle vrai : il ne suppose pas que la vérité figure parmi les candidats, et il ne devient pas certain quand l'échantillon grandit. Voir [1a], [1c].
2. [2]
   BIC, the Bayesian information criterion, is p ln N − 2 ln L. Same shape, different price: ln N instead of 2. It comes from a different question, namely an approximation to the probability the data assign to each model after integrating over its parameters, so it estimates which model generated the data rather than which will predict best. Because the penalty grows with the sample, BIC keeps demanding more evidence for each extra parameter as data accumulate, and it will settle on the true model given enough of them, if the true model is among the candidates. From 8 points on, ln N exceeds 2, so BIC is the stricter of the two for any realistic sample. See [1b], [1d].
   Le BIC, critère d'information bayésien, vaut p ln N − 2 ln L. Même forme, prix différent : ln N au lieu de 2. Il découle d'une autre question, à savoir une approximation de la probabilité que les données accordent à chaque modèle après intégration sur ses paramètres ; il estime donc quel modèle a engendré les données plutôt que lequel prédira le mieux. Comme la pénalité croît avec l'échantillon, le BIC exige toujours plus de preuves pour chaque paramètre supplémentaire à mesure que les données s'accumulent, et il finit par retenir le modèle vrai s'il figure parmi les candidats. Dès 8 points, ln N dépasse 2 : le BIC est donc le plus strict des deux pour tout échantillon réaliste. Voir [1b], [1d].
3. [3]
   Why a difference below about 2 counts as a tie. A BIC difference is roughly twice the logarithm of the Bayes factor between the two models, and on the scale Kass and Raftery [2a] set out for reading Bayes factors, that bottom band is the one they describe as not worth more than a bare mention. The same threshold is conventional for AIC. This is why select_orders gives you the whole ranking rather than only the winner: if the top few are within 2 of each other, the criterion has not chosen anything, and the simplest of them is the one to keep.
   Pourquoi un écart inférieur à environ 2 vaut une égalité. Un écart de BIC vaut approximativement deux fois le logarithme du facteur de Bayes entre les deux modèles, et sur l'échelle proposée par Kass et Raftery [2a] pour lire les facteurs de Bayes, cette première bande est celle qu'ils décrivent comme ne méritant pas plus qu'une simple mention. Le même seuil est d'usage pour l'AIC. C'est pourquoi select_orders vous renvoie tout le classement et pas seulement le gagnant : si les premiers sont à moins de 2 les uns des autres, le critère n'a rien tranché, et c'est le plus simple qu'il faut retenir.
4. [4]
   Why maximum-likelihood widths come out too narrow. The fitted trend follows the noise a little, so the residuals around it are systematically smaller than the true deviations, by a factor of about √((N − p)/N). This is the same effect that makes a plain standard deviation divide by N − 1 rather than N [3f], generalised to p fitted parameters. dof_correction=True undoes it and is the default.
   Pourquoi les largeurs par maximum de vraisemblance sont trop étroites. La tendance ajustée épouse un peu le bruit, de sorte que les résidus autour d'elle sont systématiquement plus petits que les écarts vrais, d'un facteur d'environ √((N − p)/N). C'est le même effet qui fait diviser un écart-type ordinaire par N − 1 plutôt que par N [3f], généralisé à p paramètres ajustés. dof_correction=True le corrige et constitue le comportement par défaut.
5. [5]
   t₃, and the Student-t family. nu is the number of degrees of freedom of a Student-t distribution, written t_ν, and it controls how heavy the tails are. As ν grows the distribution becomes a Gaussian: by ν = 30 the two are already hard to tell apart. As ν shrinks the tails thicken. Take t₃, the column in the table above. Its scale parameter is still called σ, but its standard deviation is √3 σ, about 1.73 times larger, and excursions that a Gaussian treats as impossible become ordinary: a 3σ departure happens 21 times more often than under a Gaussian, and a 5σ one roughly 27 000 times more often. Push further and the family stops being well behaved at all: at ν ≤ 2 the variance is infinite, and at ν = 1 it is the Cauchy distribution, which has no mean.
   t₃, et la famille de Student. nu est le nombre de degrés de liberté d'une loi de Student, notée t_ν, et il gouverne l'épaisseur des queues. Quand ν croît, la loi tend vers une gaussienne : dès ν = 30 les deux sont difficiles à distinguer. Quand ν diminue, les queues s'épaississent. Prenez t₃, la colonne du tableau ci-dessus. Son paramètre d'échelle s'appelle toujours σ, mais son écart-type vaut √3 σ, environ 1,73 fois plus, et des écarts qu'une gaussienne juge impossibles y deviennent ordinaires : un écart de 3σ y survient 21 fois plus souvent que sous une gaussienne, et un écart de 5σ environ 27 000 fois plus souvent. En poussant plus loin, la famille cesse d'être raisonnable : à ν ≤ 2 la variance est infinie, et à ν = 1 c'est la loi de Cauchy, qui n'a même pas de moyenne.

ReferencesRéférences