« Découvrir Scilab/Régression » : différence entre les versions

m
<source> -> <syntaxhighlight> (phab:T237267)
m (typo)
m (<source> -> <syntaxhighlight> (phab:T237267))
: Y = (''y''<sub>1</sub> ''y''<sub>2</sub> … ''y<sub>n&thinsp;</sub>'')
Si l'on demande à Scilab de résoudre le système avec la division matricielle alors que la matrice X n'est pas carrée, il effectue automatiquement une régression linéaire des moindres carrés :
<sourcesyntaxhighlight lang="scilab">
\\ génération des données de test
x = 1:10
X = [x ; ones(x)];
A = y/X
</syntaxhighlight>
</source>
 
Si l'on veut forcer l'ordonnée à l'origine à 0, il suffit de ne pas mettre la ligne de 1 :
<sourcesyntaxhighlight lang="scilab">
A = y/x
</syntaxhighlight>
</source>
 
On peut aussi travailler avec des vecteurs colonne, et utiliser l'autre diviseur matriciel <code>A = x\y</code>.
Pour cela, on peut utiliser la fonction <code>reglin()</code>, avec la syntaxe :
 
<sourcesyntaxhighlight lang = "scilab">
// définition du nuage de points
X = […];
// régression
[a, b, sigma] = reglin(X, Y);
</syntaxhighlight>
</source>
 
{{note
On peut aussi utiliser la commande <code>regress()</code>. Elle ne permet pas d'obtenir l'écart type, et ne fonctionne qu'avec deux variables, mais par contre, les vecteurs peuvent être indifféremment ligne ou colonne ; ils peuvent même être de types différents (un vecteur ligne et l'autre colonne).
 
<sourcesyntaxhighlight lang = "scilab">
coefs = regress(X, Y)
</syntaxhighlight>
</source>
 
'''Exemple'''
 
: Fichier <code>fonction.sce</code>
<sourcesyntaxhighlight lang="scilab">
chdir("monrepertoire/");
 
 
deff("[y] = affine(a, b, x)", "y = a*x + b")
</syntaxhighlight>
</source>
 
: Fichier <code>generation.sce</code>
<sourcesyntaxhighlight lang="scilab">
chdir("monrepertoire/");
 
 
write("donnees.txt", [X,Y]);
</syntaxhighlight>
</source>
 
: Fichier <code>regression.sce</code>
<sourcesyntaxhighlight lang="scilab">
clf;
chdir("monrepertoire/");
+" ; b = "+string(bopt)+" pour "+string(ordonnee)...
+" ; sigma = "+string(sigma))
</syntaxhighlight>
</source>
 
Comme la forme de la fonction ƒ est connue, il est inutile de la définir. La fonction <code>reglin()</code> retourne un vecteur de trois composantes : les paramètres ''a'' et ''b'' du modèle, et l'écart type σ du résidu (qui est la racine carrée de la somme des carrés des résidus, σ = √S).
Voici par exemple un programme générant ''n'' points pour ''x'' allant de 0 à 5 (valeurs exactes), et ''y'' suivant le polynôme ''x''<sup>3</sup> - 2·''x''<sup>2</sup> + 3·''x'' - 4, bruité :
 
<sourcesyntaxhighlight lang="scilab">
// n points aléatoires approchant un polynôme de degré 3
 
 
write("polynome_bruite.txt", [x' , y']);
</syntaxhighlight>
</source>
 
Et voici un programme permettant d'effectuer la régression sur ces données :
 
<sourcesyntaxhighlight lang="scilab">
// Régression polynomiale de degré 3
 
plot(x, yth)
plot(x, y, "+")
</syntaxhighlight>
</source>
 
== Régression non-linéaire ==
 
Dans le cas de la régression non-linéaire, ƒ est une fonction quelconque. Pour effectuer la régression, on utilise la fonction <code>leastsq()</code> avec la syntaxe suivante :
<sourcesyntaxhighlight lang="scilab">
// fonction modèle
function y = f(A, x)
 
[S, Aopt] = leastsq(list(r, X, Y), Ainit)
</syntaxhighlight>
</source>
Avec cette syntaxe, la fonction <code>leastsq()</code> retourne donc un vecteur :
* le premier élément, <code>S</code>, est la valeur de la somme des carrés des écarts ;
 
Nous reprenons l'exemple de la régression linéaire, mais en utilisant cette fois-ci la fonction <code>leastsq()</code>. Lorsque l'on enregistre les données générées, on utilise :
<sourcesyntaxhighlight lang="scilab">
write("donnees.txt", [X', Y']);
</syntaxhighlight>
</source>
afin d'avoir des vecteurs colonne.
 
<sourcesyntaxhighlight lang="scilab">
clf;
chdir("monrepertoire/");
+" ; b = "+string(Aopt(2))+" pour "+string(ordonnee)...
+" ; variance = "+string(S))
</syntaxhighlight>
</source>
 
[[File:Desommation gaussienne lorentzienne.svg|vignette|upright=2|Résultat de la régression.]]
Nous travaillons ici avec un pic généré par la somme d'une fonction gaussienne et d'une fonction lorentzienne. Ce genre de profil est classiquement utilisé pour décrire la forme des pics obtenus par diffractométrie X. On fait ce que l'on appelle une désommation de pics (on parle souvent de « déconvolution », mais ce terme est abusif puisqu'il ne s'agit pas de produit de convolution).
 
<sourcesyntaxhighlight lang="scilab">
clf()
 
xstring(-demielargeur, hauteur/2,...
"variance : S = "+string(0.01*round(100*S))))
</syntaxhighlight>
</source>
 
L'algorithme par défaut de <code>leastsq()</code> est l'algorithme de pseudo-Newton. On peut aussi utiliser l'[[w:Méthode du gradient conjugué|algorithme du gradient conjugué]] en ajoutant <code>"gc"</code> aux paramètres, ou bien algorithme foncitonnant avec les fonctions non-différentiables avec <code>"nd"</code>. On peut enfin utiliser la commande <code>lsqrsolve()</code>, qui utilise l'algorithme de Levenberg-Marquardt.
1 535

modifications