Qui a créé la première table normale standard?

61

Je suis sur le point d'introduire la table normale standard dans mon cours d'introduction aux statistiques et je me suis demandé: qui a créé la première table normale standard? Comment l'ont-ils fait avant l'arrivée des ordinateurs? Je frémis en pensant à une calculatrice en force brute que mille sommes Riemann à la main.

Daniel Smolkin
la source
5
Heureux de voir quelqu'un qui souhaite avoir un enseignement historiquement informé.
mdewey

Réponses:

62

Laplace a été le premier à reconnaître le besoin de tabulation et a fourni l'approximation suivante:

G(x)=xet2dt(1)=1x12x3+134x51358x7+135716x9+

Le premier tableau moderne de la distribution normale a ensuite été construit par l'astronome français Christian Kramp dans Analyse des fractions astronomiques et terrestres (Par le citoyen Kramp, Professeur de Chymie et de Physique expérimentale à l'école centrale du Département de la Roer, 1799) . Des tableaux relatifs à la distribution normale: un bref historique Auteur (s): Herbert A. David Source: The American Statistician, Vol. 59, n ° 4 (nov. 2005), p. 309-311 :

Ambitieusement, Kramp a fourni des tableaux à huit décimales ( D) allant jusqu'à x = 1,24, 9 D à 1,50, 10 D à 1,99 et 11 D à 3,00, ainsi que les différences nécessaires pour l'interpolation. Notant les six premières dérivées de G (x), il utilise simplement un développement en série de Taylor de G (x + h) autour de G (x), avec h = .01, jusqu'au terme de h ^ 3. Cela lui permet de procéder pas à pas de x = 0 à x = h, 2h, 3h, \ dots, en multipliant h \, e ^ {- x ^ 2} par8x=1.24, 91.50, 101.99,113.00G(x),G(x+h)G(x),h=.01,h3.x=0x=h,2h,3h,,hex2

1hx+13(2x21)h216(2x33x)h3.
Ainsi, àx=0 ce produit est réduit à
.01(113×.0001)=.00999967,
sorte que pourG(.01)=.88622692.00999967=.87622725.


entrez la description de l'image ici

entrez la description de l'image ici

entrez la description de l'image ici

Mais ... quelle précision pourrait-il avoir? OK, prenons 2.97 comme exemple:

entrez la description de l'image ici

Incroyable!

Passons à l'expression moderne (normalisée) du pdf gaussien:

Le pdf de N(0,1) est:

fX(X=x)=12πex22=12πe(x2)2=12πe(z)2

z=x2 . Et par conséquent, x=z×2 .

Alors, allons à R et cherchons le ... OK, pas si vite. Nous devons d’abord nous rappeler que s’il existe une constante qui multiplie l’exposant dans une fonction exponentielle , l’intégrale sera divisée par cet exposant: . Puisque nous essayons de reproduire les résultats dans les anciennes tables, nous multiplions en fait la valeur de par , ce qui devra apparaître dans le dénominateur.PZ(Z>z=2.97)eax1/ax2

De plus, Christian Kramp n’ayant pas normalisé, nous devons donc corriger les résultats donnés par R en conséquence, en multipliant par . La correction finale ressemblera à ceci:2π

2π2P(X>x)=πP(X>x)

Dans le cas ci-dessus, et . Maintenant passons à R:z=2.97x=z×2=4.200214

(R = sqrt(pi) * pnorm(x, lower.tail = F))
[1] 0.00002363235e-05

Fantastique!

Allons au sommet de la table pour le plaisir, disons ...0.06

z = 0.06
(x = z * sqrt(2))

(R = sqrt(pi) * pnorm(x, lower.tail = F))
[1] 0.8262988

Que dit Kramp? .0.82629882

Si proche ...


La chose est ... à quelle distance, exactement? Après tous les votes positifs reçus, je ne pouvais pas laisser la réponse réelle en suspens. Le problème était que toutes les applications de reconnaissance optique de caractères (OCR) que j'avais essayées étaient incroyablement inactives - ce qui n'a rien d'étonnant si vous avez jeté un coup d'œil à l'original. J'ai donc appris à apprécier Christian Kramp pour la ténacité de son travail, car j'ai personnellement saisi chaque chiffre dans la première colonne de son tableau Première .

Après une aide précieuse de @Glen_b, il peut maintenant être très précis et il est prêt à copier et coller sur la console R dans ce lien GitHub .

Voici une analyse de la précision de ses calculs. Préparez vous...

  1. Différence cumulative absolue entre les valeurs [R] et l'approximation de Kramp:

0.000001200764 - en calculs, il a réussi à accumuler une erreur d'environ millionième!3011

  1. Erreur absolue moyenne (MAE) , oumean(abs(difference))avecdifference = R - kramp:

0.000000003989249 - il a réussi à commettre une erreur scandaleusement ridicule au milliardième en moyenne!3

Sur l'entrée dans laquelle ses calculs étaient les plus divergents par rapport à [R], la première décimale différente se situait à la huitième position (cent millionième). En moyenne (médiane), sa première "erreur" était dans le dixième chiffre décimal (dixième milliardième!). Et, bien qu'il ne soit absolument pas d'accord avec [R] dans aucun cas, l'entrée la plus proche ne divergent que jusqu'à la 13 entrée numérique.

  1. Différence relative moyenne ou mean(abs(R - kramp)) / mean(R)(identique à all.equal(R[,2], kramp[,2], tolerance = 0)):

0.00000002380406

  1. Erreur quadratique moyenne (RMSE) ou écart (donne plus de poids aux erreurs importantes), calculé comme suitsqrt(mean(difference^2)):

0.000000007283493


Si vous trouvez une photo ou un portrait de Chistian Kramp, éditez ce message et placez-le ici.

Antoni Parellada
la source
4
C'est bien d'avoir les deux références différentes, et je pense que les détails supplémentaires (comme l'extension explicite donnée par Laplace pour la queue supérieure) sont bons.
Glen_b
1
C'est encore mieux avec le dernier montage, mais je ne peux pas voter deux fois, c'est excellent. Notez que l'article de David explique pourquoi la table de Kramp ne précise pas tous les chiffres indiqués (une très petite erreur a été commise dans la première étape) - mais cela reste plus que suffisant pour la plupart des applications statistiques
Glen_b 5/05/2016
2
@ OlivierGrégoire Merci d'avoir signalé mon chiffre décimal mal orthographié. Il est maintenant corrigé. J'ai grandi à une époque où le français était une nécessité et ne signifiait en aucun cas un manque de respect vis-à-vis de mon utilisation bizarre de la langue (il y a une référence dedans, mais peu importe), que j'ai inversée. En ce qui concerne "citoyen Kramp" - une tentative de mise en évidence des formes historiques d'introduction dans le journal.
Antoni Parellada
1
Hé, désolé que tu aies senti que c'était un commentaire bashing. Je ne faisais que pointer du doigt, je ne vous dis absolument pas que vous manquiez de respect. Vous pouvez punir ou exagérer (ou même faire une référence), bien sûr. Mais en tant que francophone, je ne l’ai pas compris (c’est ce que j’ai essayé de transmettre, du moins). "Le citoyen Kramp" n'avait pas de problème: je venais de copier et de mettre des guillemets, car ce n'était pas l'anglais. Désolé si vous sentiez que c'était un commentaire bashing, ce n'est pas. Mon utilisation de l'anglais fait également défaut. ^^ Votre comparaison a été bien faite!
Olivier Grégoire
1
@ P.Windridge Désolé ... j'ai réalisé que j'avais un tas de liens hypertextes cassés ...
Antoni Parellada
32

Selon HA, David [1] Laplace a reconnu la nécessité de disposer de tables de distribution normale "dès 1783" et la première table normale a été produite par Kramp en 1799.

Laplace a suggéré deux approximations en série, une pour l'intégrale de à de (qui est proportionnelle à une distribution normale avec variance ) et une pour la queue supérieure.x e - t 2 10xet212

Cependant, Kramp n'a pas utilisé ces séries de Laplace, car il y avait un intervalle dans les intervalles pour lesquelles elles pourraient être utilement appliquées.

En fait, il part de 0 pour l'intégrale de la zone de la queue, puis applique une expansion de Taylor sur la dernière intégrale calculée. En calculant de nouvelles valeurs dans le tableau, il décale le de son expansion de Taylor de (où est l'intégrale donnant la partie supérieure de la queue).G ( x + h ) GxG(x+h)G

Pour être précis, citant les quelques phrases pertinentes:

il utilise simplement un développement en série de Taylor de autour de , avec , jusqu'au terme de . Cela lui permet de procéder pas à pas de à , en multipliant parAinsi, à ce produit est réduit à sorte que, à . On peut montrer que le terme suivant sur la gauche de (4) est , de sorte que son omission est justifiée.G(x+h)G(x)h=.01h3x=0x=h,2h,3h,...hex2x=00,01(1-1

1hx+13(2x21)h216(2x33x)h3.
x=0G ( .01 ) = .88622692 - .00999967 = .87622725 10 - 9
.01(113×.0001)=.00999967,(4)
G(.01)=.88622692.00999967=.87622725109

David indique que les tables ont été largement utilisées.

Ainsi, plutôt que des milliers de sommes Riemann, il s’agissait de centaines d’agrandissements de Taylor.


Sur une note plus petite, avec un pincement (coincé avec seulement une calculatrice et quelques valeurs mémorisées du tableau normal), j'ai assez bien appliqué la règle de Simpson (et les règles associées pour l'intégration numérique) pour obtenir une bonne approximation à d'autres valeurs; ce n'est pas tout ce que fastidieux de produire une table abrégée * à quelques chiffres de précision. [Produire des tableaux d'échelle et de précision de Kramp serait une tâche assez lourde, même avec une méthode plus intelligente, comme il l'a fait.]

* Par tableau abrégé, je veux dire un tableau dans lequel vous pouvez vous échapper avec une interpolation entre les valeurs tabulées sans perdre trop de précision. Si vous voulez seulement dire une précision d'environ 3 chiffres, vous n'avez vraiment pas besoin de calculer autant de valeurs. J'ai utilisé efficacement l'interpolation polynomiale (plus précisément, les techniques de différences finies appliquées), ce qui permet d'obtenir un tableau avec moins de valeurs que l'interpolation linéaire - même si l'effort requis à l'étape d'interpolation est plus important - et j'ai également effectué une interpolation avec une transformation logit, qui rend l’interpolation linéaire considérablement plus efficace, mais n’est utile que si vous avez une bonne calculatrice).

[1] Herbert A. David (2005),
"Tableaux relatifs à la distribution normale: une courte histoire",
The American Statistician , vol. 59, n ° 4 (nov.), P. 309-311

[2] Kramp (1799),
Analyse des fractions astronomiques et terrestres,
Leipzig: Schwikkert.

Glen_b
la source
0

Question intéressante! Je pense que la première idée n’est pas venue de l’intégration de formules complexes; plutôt le résultat de l'application des asymptotiques en combinatoire. La méthode du stylo et du papier peut prendre plusieurs semaines; pas si difficile pour Karl Gauss par rapport au calcul du gâteau pour ses prédécesseurs. Je pense que l'idée de Gauss était courageuse; le calcul était facile pour lui.

Exemple de création d'une table z standard à partir de zéro-
1. Prenez une population de n nombres (disons que n vaut 20) et dressez une liste de tous les échantillons possibles de taille r (disons que r est égale à 5).
2. calculer les moyennes d'échantillon. Vous obtenez des moyennes d'échantillon nCr (ici, 20c5 = 15504 signifie).
3. Leur moyenne est la même que la moyenne de la population. Trouver le stdev des moyennes d'échantillon.
4. Trouvez les scores z des moyennes d'échantillon en utilisant ces moyennes pop et stdev des moyennes d'échantillon.
5. Triez z par ordre croissant et trouvez la probabilité que z se situe dans une plage de vos valeurs nCr z.
6. Comparez les valeurs avec des tables normales. Plus petit n est bon pour les calculs manuels. Plus grand n produira des approximations plus proches des valeurs normales de la table.

Le code suivant est en r:

n <- 20  
r <- 5  

p <- sample(1:40,n)  # Don't be misled!! Here, 'sample' is an r function  
                     used to produce n random numbers between 1 and 40.  
                     You can take any 20 numbers, possibly all different.  

c <- combn(p, r)     # all the nCr samples listed  
cmean <- array(0)  

for(i in 1:choose(n,r)) {  
    cmean[i] <- mean(c[,i])  
                }  

z <- array(0)  
for(i in 1:choose(n,r)) {  
    z[i] <- (cmean[i]-mean(c))/sd(cmean)  
                }  

ascend <- sort(z, decreasing = FALSE)  

Probabilité que z tombe entre 0 et une valeur positive q inférieure; comparer avec un tableau connu. Manipulez q ci-dessous entre 0 et 3,5 pour comparer.

q <- 1  
probability <- (length(ascend[ascend<q])-length(ascend[ascend<0]))/choose(n,r)   
probability   # For example, if you use n=30 and r=5, then for q=1, you  
              will get probability is 0.3413; for q=2, prob is 0.4773
Md Towhidul Islam
la source
3
Je ne vois pas comment l’échantillonnage peut être utilisé de cette manière pour générer les tableaux. Je pense que l'OP voulait juste savoir qui était la première personne
Michael Chernick
Merci pour votre précieux commentaire, Michael Chernick. 1) Le PO écrit: "Comment ont-ils procédé avant l’arrivée des ordinateurs? Je frémis en pensant à une calculatrice en force brute mille dollars Riemann à la main." J'ai essayé de répondre à cette partie. 2) Le terme «échantillon» n’est pas un échantillon en soi, c’est une fonction r qui consiste à produire une liste de nombres aléatoires. Nous pouvons simplement prendre 20 numéros à la place. Voir le lien de soutien ici stackoverflow.com/questions/17773080/…
Md Towhidul Islam