Comment trouver le mode d'une fonction de densité de probabilité?

Inspiré par mon autre question , j'aimerais savoir comment trouver le mode d'une fonction de densité de probabilité (PDF) d'une fonction ?$f(x)$

Existe-t-il une procédure de "livre de cuisine" pour cela? Apparemment, cette tâche est beaucoup plus difficile qu'il n'y paraît au début.

Dire "le mode" implique que la distribution en a un et un seul. En général, une distribution peut avoir de nombreux modes, ou (sans doute) aucun.

S'il y a plus d'un mode, vous devez spécifier si vous voulez tous ou juste le mode global (s'il y en a exactement un).

En supposant que nous nous limitions aux distributions unimodales *, afin que nous puissions parler du "mode", elles se retrouvent de la même manière que la recherche de maxima de fonctions plus généralement.

* notez que la page dit " comme le terme" mode "a plusieurs significations, ainsi que le terme" unimodal " " et propose plusieurs définitions de mode - qui peuvent changer ce qui, exactement, compte comme mode, qu'il y ait 0 1 ou plus - et modifie également la stratégie pour les identifier. Notez en particulier à quel point le libellé "plus général" de ce qu'est l'unimodalité est général dans le premier paragraphe " unimodalité signifie qu'il n'y a qu'une seule valeur la plus élevée, en quelque sorte définie "

Une définition proposée sur cette page est:

Un mode de distribution de probabilité continue est une valeur à laquelle la fonction de densité de probabilité (pdf) atteint sa valeur maximale

Donc, étant donné une définition spécifique du mode, vous le trouvez comme vous trouveriez cette définition particulière de "valeur la plus élevée" lorsque vous traitez les fonctions de manière plus générale (en supposant que la distribution est unimodale selon cette définition).

Il existe une variété de stratégies en mathématiques pour identifier de telles choses, selon les circonstances. Voir la section «Trouver les maxima et les minima fonctionnels» de la page Wikipedia sur les maxima et les minima qui donne une brève discussion.

Par exemple, si les choses sont suffisamment belles - disons que nous avons affaire à une variable aléatoire continue, où la fonction de densité a une dérivée première continue - vous pouvez continuer en essayant de trouver où la dérivée de la fonction de densité est nulle et en vérifiant de quel type de point critique il s'agit (maximum, minimum, point d'inflexion horizontal). S'il y a exactement un tel point qui est un maximum local, ce devrait être le mode d'une distribution unimodale.

Cependant, en général, les choses sont plus compliquées (par exemple, le mode peut ne pas être un point critique), et les stratégies plus larges pour trouver des maxima de fonctions entrent en jeu.

Parfois, trouver algébriquement où les dérivés sont nuls peut être difficile ou au moins lourd, mais il peut toujours être possible d'identifier les maxima par d'autres moyens. Par exemple, il se peut que l'on puisse invoquer des considérations de symétrie pour identifier le mode d'une distribution unimodale. Ou on pourrait invoquer une certaine forme d'algorithme numérique sur un ordinateur, pour trouver un mode numériquement.

Voici quelques cas qui illustrent des choses typiques que vous devez vérifier - même lorsque la fonction est unimodale et au moins continue par morceaux.

Ainsi, par exemple, nous devons vérifier les points de terminaison (diagramme central), les points où le dérivé change de signe (mais peut-être pas zéro; premier diagramme) et les points de discontinuité (troisième diagramme).

Dans certains cas, les choses peuvent ne pas être aussi soignées que ces trois-là; vous devez essayer de comprendre les caractéristiques de la fonction particulière à laquelle vous avez affaire.

Je n'ai pas abordé le cas multivarié, où même lorsque les fonctions sont assez "agréables", trouver des maxima locaux peut être considérablement plus complexe (par exemple, les méthodes numériques pour le faire peuvent échouer dans un sens pratique, même lorsqu'elles doivent logiquement réussir) finalement).

Cette réponse se concentre entièrement sur l'estimation de mode à partir d'un échantillon, en mettant l'accent sur une méthode particulière. S'il y a un sens fort dans lequel vous connaissez déjà la densité, analytiquement ou numériquement, alors la réponse préférée est, en bref, de rechercher directement le maximum unique ou les maxima multiples, comme dans la réponse de @Glen_b.

Les "modes de demi-échantillon" peuvent être calculés en utilisant une sélection récursive du demi-échantillon de la plus courte longueur. Bien qu'elle ait des racines plus longues, une excellente présentation de cette idée a été donnée par Bickel et Frühwirth (2006).

L'idée d'estimer le mode comme le milieu de l'intervalle le plus court contenant un nombre fixe d'observations remonte au moins à Dalenius (1965). Voir aussi Robertson et Cryer (1974), Bickel (2002) et Bickel et Frühwirth (2006) sur d'autres estimateurs du mode.

Les statistiques d'ordre d'un échantillon de valeurs de x sont définies par x ( 1 ) ≤ x ( 2 ) ≤ ⋯ ≤ x ( n - 1 ) ≤ x ( n ) .$n$$x$$x_{(1)} \le x_{(2)} \le \cdots \le x_{(n-1)} \le x_{(n)}$

Le mode demi-échantillon est ici défini à l'aide de deux règles.

Règle 1. Si , le mode demi-échantillon est x ( 1 ) . Si n = 2 , le mode demi-échantillon est ( x ( 1 ) + x ( 2 ) ) / 2 . Si n = 3 , le mode demi-échantillon est ( x ( 1 ) + x ( 2 ) ) / 2 si x ( 1 ) et x ( $n = 1$$x_{(1)}$$n = 2$$(x_{(1)} + x_{(2)}) / 2$$n = 3$$(x_{(1)} + x_{(2)}) / 2$$x_{(1)}$$x_{(2)}$$x_{(2)}$$x_{(3)}$$(x_{(2)} + x_{(3)}) / 2$$x_{(2)}$

$n \ge 4$$3$$h_1 = \lfloor n / 2\rfloor$$k$$k + h_1$$x_{(k + h_1)} - x_{(k)}$$k = 1, \cdots, n - h_1$$h_1 + 1$$h_2 = \lfloor h_1 / 2\rfloor$

$x_{(k)}, \cdots, x_{(k + h)}$$h = \lfloor n / 2 \rfloor$ est l'estimateur le moins médian des carrés (LMS) de l'emplacement pour x . Voir Rousseeuw (1984) et Rousseeuw et Leroy (1987) pour les applications du LMS et les idées associées à la régression et à d'autres problèmes. Notez que ce point médian du LMS est également appelé le shorth dans certains ouvrages plus récents (par exemple Maronna, Martin et Yohai 2006, p.48). De plus, la moitié la plus courte elle-même est parfois appelée shorth, comme l'indique le titre de Grübel (1988). Pour une implémentation Stata et plus de détails, voir à partir de SSC.$(x_k + x_{(k + h)}) / 2$$x$shorth

Certains commentaires généraux suivent les avantages et les inconvénients des modes de demi-échantillon, du point de vue des analystes de données pratiques autant que des statisticiens mathématiques ou théoriques. Quel que soit le projet, il sera toujours judicieux de comparer les résultats avec des mesures sommaires standard (par exemple des médianes ou des moyennes, y compris des moyennes géométriques et harmoniques) et de relier les résultats à des graphiques de distributions. De plus, si vous vous intéressez à l'existence ou à l'étendue de la bimodalité ou de la multimodalité, il sera préférable de regarder directement des estimations convenablement lissées de la fonction de densité.

Estimation de mode En résumant où les données sont les plus denses, le mode demi-échantillon ajoute un estimateur automatique du mode à la boîte à outils. Des estimations plus traditionnelles du mode basées sur l'identification de pics sur des histogrammes ou même des graphiques de densité de noyau sont sensibles aux décisions sur l'origine ou la largeur du bac ou le type et la demi-largeur du noyau et sont plus difficiles à automatiser de toute façon. Lorsqu'il est appliqué à des distributions unimodales et approximativement symétriques, le mode demi-échantillon sera proche de la moyenne et de la médiane, mais plus résistant que la moyenne aux valeurs aberrantes dans les deux queues. Lorsqu'il est appliqué à des distributions unimodales et asymétriques, le mode demi-échantillon sera généralement beaucoup plus proche du mode identifié par d'autres méthodes que la moyenne ou la médiane.

Simplicité L'idée du mode demi-échantillon est assez simple et facile à expliquer aux étudiants et chercheurs qui ne se considèrent pas comme des statisticiens.

Interprétation graphique Le mode demi-échantillon peut facilement être lié à des affichages standard de distributions telles que les diagrammes de densité de noyau, les diagrammes de distribution cumulative et quantile, les histogrammes et les diagrammes tige-feuille.

En même temps, notez que

Inutile pour toutes les distributions Lorsqu'il est appliqué à des distributions qui sont approximativement en forme de J, le mode demi-échantillon approximera le minimum des données. Lorsqu'il est appliqué à des distributions qui sont approximativement en forme de U, le mode demi-échantillon sera dans la moitié de la distribution qui a la densité moyenne la plus élevée. Aucun des deux comportements ne semble particulièrement intéressant ou utile, mais de même, il n'y a guère de demande de résumés de type monomode pour les distributions en forme de J ou en U. Pour les formes en U, la bimodalité rend l'idée d'un moot mode unique, sinon invalide.

Liens La moitié la plus courte peut ne pas être définie de manière unique. Même avec des données mesurées, l'arrondi des valeurs déclarées peut fréquemment donner lieu à des liens. Que faire des deux moitiés les plus courtes ou plus a été peu discuté dans la littérature. Notez que les moitiés liées peuvent se chevaucher ou être disjointes.

hsmode$t$$t$$\lceil t/ 2\rceil$

$-9, -4, -1 , 0, -1, 4, 9$$-0.5$$0$$1 + \lfloor n / 2\rfloor$$n$$n$, ce qui est difficile à réaliser étant donné d'autres desiderata, notamment que la longueur de la fenêtre ne devrait jamais diminuer avec la taille de l'échantillon. Nous préférons croire qu'il s'agit d'un problème mineur avec des ensembles de données de taille raisonnable.

$1 + \lfloor n / 2\rfloor$$n$$n$$n = 1,$$n = 2$$\lceil n / 2\rceil$

$1.6, 3.11, 3.95, 4.2, 4.2, 4.62, 4.62, 4.62, 4.7, 4.87, 5.04, 5.29, 5.3, 5.38, 5.38, 5.38, 5.54, 5.54, 5.63, 5.71, 6.13, 6.38, 6.38, 6.67, 6.69, 6.97, 7.22, 7.72, 7.98, 7.98, 8.74, 8.99, 9.27, 9.74, 10.66.$hsmode$5.00, 5.02, 5.04$

Andrews, DF, PJ Bickel, FR Hampel, PJ Huber, WH Rogers et JW Tukey. 1972. Estimations robustes de l'emplacement: relevé et progrès. Princeton, NJ: Princeton University Press.

Bickel, DR 2002. Estimateurs robustes du mode et de l'asymétrie des données continues. Statistiques computationnelles et analyse des données 39: 153-163.

Bickel, DR et R. Frühwirth. 2006. Sur un estimateur rapide et robuste du mode: comparaisons avec d'autres estimateurs avec applications. Statistiques computationnelles et analyse des données 50: 3500-3530.

Dalenius, T. 1965. Le mode - Un paramètre statistique négligé. Journal, Royal Statistical Society A 128: 110-117.

Grübel, R. 1988. La longueur du shorth. Annals of Statistics 16: 619-628.

Hampel, FR 1975. Au-delà des paramètres de localisation: concepts et méthodes robustes. Bulletin, Institut international de statistique 46: 375-382.

Maronna, RA, RD Martin et VJ Yohai. 2006. Statistiques robustes: théorie et méthodes . Chichester: John Wiley.

Robertson, T. et JD Cryer. 1974. Une procédure itérative pour estimer le mode. Journal, American Statistical Association 69: 1012-1016.

Rousseeuw, PJ 1984. La médiane de régression des moindres carrés. Journal, American Statistical Association 79: 871-880.

Rousseeuw, PJ et AM Leroy. 1987. Régression robuste et détection des valeurs aberrantes . New York: John Wiley.

Ce compte est basé sur la documentation de

Cox, NJ 2007. HSMODE: module Stata pour calculer les modes de demi-échantillon, http://EconPapers.repec.org/RePEc:boc:bocode:s456818 .

Voir également le site Web de David R. Bickel ici pour des informations sur les implémentations dans d'autres logiciels.

Si vous avez des échantillons de la distribution dans un vecteur "x", je ferais:

 mymode <- function(x){
   d<-density(x)
   return(d$x[which(d$y==max(d$y)[1])])
 }

Vous devez régler la fonction de densité afin qu'elle soit suffisamment lisse sur le dessus ;-).

Si vous n'avez que la densité de la distribution, j'utiliserais un optimiseur pour trouver le mode (REML, LBFGS, simplex, etc.) ...

 fx <- function(x) {some density equation}
 mode <- optim(inits,fx)

Ou utilisez un échantillonneur Monte-Carlo pour obtenir des échantillons de la distribution (package rstan) et utilisez la procédure ci-dessus. (Quoi qu'il en soit, le package Stan comme une fonction "d'optimisation" pour obtenir le mode d'une distribution).