Découvrir Scilab/Statistiques
4. Statistiques
Pour l'analyse statistique de données, on peut se reporter au document
Tendance centrale.
Considérons une matrice X. On peut calculer les moyennes suivantes :
On peut calculer la moyenne géométrique avec codice_4, et la moyenne harmonique avec codice_5, avec la même syntaxe.
De la même manière, on peut calculer les médianes avec codice_6, codice_7 et codice_8.
La commande codice_9 permet de calculer la moyenne en éliminant les "p" pourcents extrêmes (les "p"/2 % valeurs les plus élevées, et les "p"/2 % valeurs les plus faibles). On peut également faire le calcul par ligne ou par colonne — codice_10 ou codice_11. notons que l'on a :
Dispersion.
La commande
strange(X)
calcule l'écart "(range)" entre la valeur maximale et la valeur minimale de la matrice. On peut calculer les valeurs pour chaque colonne ou chaque ligne avec codice_14 et codice_15. 
La commande
mad(X)
calcule la déviation absolue moyenne "(mean absolute deviation)"
on peut de même la calculer pour chaque colonne ou ligne.
La commande
iqr(X)
calcule l'écart interquartile "(interquartile range)".
La commande codice_16 calcule l'écart type "(standard deviation)" et la commande codice_17 la variance des éléments de la matrice (entière, par ligne ou par colonne). Par défaut, la variance est sans biais : normalisée par "m" × "n" – 1 pour une matrice de "m" × "n" éléments, par "m" – 1 par colonne et par "n" – 1 par ligne. Si l'on calcule la variance par ligne ou par colonne, on peut ajouter l'option codice_18 pour calculer le moment d'ordre 2 (normalisation par "m" ou "n"), par exemple
mom2 = variance(X, "r", 1) // pour une matrice non-vecteur
mom2 = variance(X, "*", 1) // pour un vecteur
L'écart type est toujours sans biais. Si l'on veut l'estimateur normalisé par "m" ou "n", il faut utiliser
estim_ecart_type = sqrt(variance(X, "*", 1)) // pour une matrice non-vecteur
Normalisation.
La commande codice_19 permet de centrer la matrice (équivalent de codice_20). La commande codice_21 donne la matrice centrée réduite (équivalent de codice_22).
Fréquences et quartiles.
La commande codice_23 calcule les fréquences d'apparition des valeurs dans la matrice X. Le résultat est classé par défaut par ordre décroissant des valeurs ; pour avoir un classement croissant, on ajoute le paramètre codice_24 "()" :
tabul(X, "i")
La commande codice_25 calcule les quartiles et les renvoie dans un vecteur de dimension 3. On peut calculer les quartiles pour chaque colonne ou chaque ligne, avec codice_26 et codice_27.
On peut aussi travailler avec des valeurs et leur fréquence. Si la matrice X contient les valeurs et la matrice F leur fréquence — F("i", "j") est la fréquence de X("i", "j") —, alors
Notons que F n'est pas nécessairement normalisé à 1 ou à 100 ; il peut par exemple s'agir d'un nombre d'occurrences.
Exemple
X = int(10*rand(100,1, "normal")); // génère 100 valeurs entières avec une loi normale centrée d'écart type 10
F = tabul(X); // calcule les fréquences d'apparition
x = F(:, 1); // extraction de la liste des valeurs
f = F(:, 2); // extraction de la liste des fréquences
quartiles = quart(X); // calcule les quartiles de la population
moyenne = meanf(x, f); // calcul de la moyenne pondérée par les fréquences
Considérons maintenant deux vecteurs X et Y de même dimension, et la matrice F des fréquences de X et Y : F("i", "j") est la fréquence de X("i") & Y("j"). Alors :
Exemple
n = 10; // nombre de valeurs
X = 1:n; // abscisses
Y = 2*X + 5 + 3*rand(X); // ordonnées : droite "bruitée"
plot(X, Y); // tracé des données
F = eye(n, n); // couplage (x, y) : chaque couple est présent une fois, (X(i) & Y(i)) = 1, (X(i) & Y(j)) = 0
covariance1 = covar(X, Y, F);
covariance2 = 1/n*sum(X.*Y) - 1/n^2*sum(X)*sum(Y); // calcul "à la main" pour comparer
correlation1 = correl(X, Y, F);
correlation2 = covariance2/sqrt(variance(X, "*", 1))/sqrt(variance(Y, "*", 1)); // calcul "à la main"
Fonctions de distribution et de répartition de lois statistiques.
Par défaut, Scilab propose une loi de distribution, pour la loi binomiale, et onze lois de répartition (probabilités cumulatives).
Pour la loi binomiale, la fonction de répartition pour "n" épreuves (nombre entier) avec une probabilité "p" (nombre réel entre 0 et 1) est donnée par
f = binomial(p, n)
La variable codice_33 est un vecteur de "n" + 1 composantes, avec
Notez que les coefficients binomiaux C"k""n" s'obtiennent eux avec la commande codice_34.
La fonction de répartition s'obtient avec codice_35 :
P = = cdfbin("PQ", k, n, p, 1-p)
donne la probabilité P d'avoir entre 0 et "k" succès pour "n" tirages avec une probabilité de réussite "p". On peut d'ailleurs le comparer avec
f = binomial(p, n)
c = sum(f(1:k+1))
De manière générale, les fonctions de répartition commencent par codice_36 "(cumulative distribution function)", par exemple :
Notons que pour chaque cas, on peut aussi avoir la valeur complémentaire Q = 1 – P avec la syntaxe
 [P, Q] = cdf…("PQ", …)
Exemple
Ci-dessous, nous traçons la fonction de répartition de la loi normale centrée réduite (μ = 0, σ = 1). Notons que pour des calculs vectoriels, l'espérance et l'écart type sont des vecteurs de même dimension que X.
X = -3:0.1:3;
mu = zeros(X);
sigma = ones(X);
[P, Q] = cdfnor("PQ", X, mu, sigma);
plot(X', [P', Q'])
Les mêmes fonctions permettent de rechercher les quantiles ou fractiles :
Exemple
Ci dessous, nous déterminons le quantile pour un risque α bilatéral de 5 %, c'est-à-dire que l'on rejette 0,25 % des cas les plus hauts et les plus bas, pour une loi normale centrée réduite puis pour une loi t de Student à 30 degrés de liberté
-->x = cdfnor("X", 0, 1, 0.975, 0.025)
 x =
 1.959964
-->t = cdft("T", 30, 0.975, 0.025)
 t =
 2.0422725
Enfin, ces fonctions permettent de déterminer les paramètres des loi pour atteindre un quantile visé :
Générateurs pseudo-aléatoires.
Les générateurs de nombres aléatoires permettent de tester des programmes ou bien de faire des simulations par la méthode de Monte Carlo.
Scilab possède une fonction codice_53 qui par défaut génère un nombre entre 0 et 1 avec une loi de distribution uniforme. Par exemple, pour simuler un jet de dé à six faces, on peut utiliser :
de = 1 + int(6*rand())
La fonction peut générer une matrice "m" × "n" avec les syntaxes
rand(m, n)
rand(X)
où X est une matrice "m" × "n". Notons que si l'on veut générer "n" valeurs, il faut bien taper codice_54 ou codice_55.
On peut lui faire générer une valeur selon une loi normale centrée réduite avec la clef codice_56 :
rand(m, n, "normal")
Exemple
Supposons par exemple que l'on a une poutre faite dans un acier dont la résistance intrinsèque, la limite d'élasticité Re, est garantie au minimum à . En réalité, Re suit une loi normale d'espérance μRe = et d'écart type σRe = . Dans les conditions d'utilisation, l'acier de cette poutre doit résister à une contrainte S "(stress)" de ; en fait, les conditions de chargement sont variables et suivent une loi normale de paramètres μS = et σS = .
On veut déterminer statistiquement dans combien de cas la contrainte S est supérieure à la résistance Re.
n = 1e6; // nombre d'événements
R = 269 + 17*rand(n, 1, "normal"); // résistance
S = 214 + 14*rand(n, 1, "normal"); // contrainte
Mbool = S>R; // comparaison des valeurs (matrice de booléens)
moyenne = sum(Mbool)/n; // les "vrai" comptent pour 1 et les "faux" pour 0
disp(string(100*moyenne)+" %") // affichage du résultat
La commande codice_57 "(get random number)" permet d'avoir des nombres aléatoires selon de nombreuses lois, par exemple :
Dans l'exemple précédent, on aurait donc :
R = grand(n, 1, "nor", 369, 17); // résistance
S = grand(n, 1, "nor", 314, 14); // contrainte
Autres lois statistiques.
Des modules complémentaires Atoms permettent d'avoir accès à d'autres loi.
Par exemple le module CASCI donne accès à de nombreuses fonctions de distribution et de répartition, notamment des lois log-normale (codice_62, codice_63), exponentielle (codice_64, codice_65) ou de Weibull (codice_66, codice_67). Ainsi :
p = pdfnormal(x, mu, sigma); // densité d'une loi normale en x
y = cdflognormal(x, lambda); // répartition d'une loi log-normale d'espérance mu et d'écart type sigma
y = cdfexponential(x, mu, sigma); // répartition d'une loi exponentielle de paramètre lambda
y = cdfweibull(x, k, lambda, theta); // répartition d'une loi de Weibull
 // de paramètre de forme k, de paramètre d'échelle lambda et de paramètre de position thêta
Le module fournit également le calcul des quantiles avec les fonctions codice_68 "()", et des générateurs de nombres aléatoires avec les fonctions codice_69. Par exemple :
x = idfweibull(y, k, lambda, theta); // quantile y d'une loi de Weibull
 // de paramètre de forme k, de paramètre d'échelle lambda et de paramètre de position thêta
X = rndweibull(n, k, lambda, theta); // génère n nombres aléatoires selon la loi de Weibull
Module Stixbox.
Le module Stixbox définit également un certain nombre de fonctions statistiques utiles. Citons par exemple la fonction codice_70 qui permet de calculer les quantiles expérimentaux et la fonction codice_71 qui permet de tracer un diagramme quantile-quantile (comme une droite de Henry). Par exemple :
n = 100;
X = rand(1, n, "norm"); // tirage aléatoire de n valeurs
// avec une loi normale centrée réduite
p = 0.1:0.1:0.9;
q = 1 - p;
qexp = quantile(X, p); // déciles expérimentaux
for i = 1:size(p, "*") // déciles théorique pour une loi normale centrée réduite
 qth(i) = cdfnor("X", 0, 1, p(i), q(i));
end
scf(0);
clf();
qqplot(qth', qexp); // diagramme quantile-quantile
xtitle("diagramme quantiles-quantiles", "quantiles théoriques"...
 "quantiles expérimentaux");
Notes.
"Calcul numérique" < ↑ > "Interpolation, extrapolation et lissage"