[Matlab] Régression non linéaire

Tout ce qui concerne notamment les outils de calcul numérique, de calcul formel ou de géométrie.
[participation réservée aux membres inscrits]
Règles du forum
Merci de soigner la rédaction de vos messages et de consulter ce sujet avant de poster. Pensez également à utiliser la fonction recherche du forum.
nirosis
Administrateur
Administrateur
Messages : 1803
Inscription : samedi 28 mai 2005, 14:48
Localisation : Orsay, France

[Matlab] Régression non linéaire

Message non lu par nirosis »

Mon problème est le suivant :

J'ai un nuage de point dans l'espace en 3 dimensions et je souhaite trouver la surface de regression (de forme polynomiale par exemple) et la visualiser avec Matlab.
J'ai du mal à comprendre les modèles dans Matlab et comment utiliser les fonctions.

J'arrive à faire des regressions en 2-D (par exemple, trouver le polynôme qui approxime au mieux mon nuage de point) mais en plusieurs dimensions je coince.

Quelqu'un a déjà fait ça ?

Merci d'avance
Kuja

Message non lu par Kuja »

Bonjour,

en régression, le modélè qu'on considère est
<center>$Y=X\beta+\epsilon$</center>
où $Y={}^t(Y_1,...,Y_n)$ est le vecteur des $n$ observations, $X$ la matrice de taille $n\times p$ des $p$ régresseurs et $\epsilon$ le vecteur $n\times 1$ des erreurs. Dans la matrice des régresseurs $X$, on peut mettre ce qu'on veut. Quelques exemples :

1) Si on suppose que $Y$ est fonction d'une seule variable $x$ dont on a $n$ réalisations $(x_1,...,x_n)$ , alors en posant

<center>$X=\left[ \begin{array}{ccccc}
1 & x_1 & x_1^2 & ... & x_1^q\\
1 & x_2 & x_2^2 & ... & x_2^q\\
\vdots & \vdots & \vdots & \vdots & \vdots\\
1 & x_n & x_n^2 & ... & x_n^q\\
\end{array}\right]$</center>
on fitte un polynôme d'ordre $q$.

2) Si on suppose que $Y$ est fonction de deux variables $u$ et $v$ dont on a $n$ réalisations $(u_1,...,u_n)$ et $(v_1,...,v_n)$ , alors en posant

<center>$X=\left[ \begin{array}{cccccc}
1 & u_1 & v_1 & u_1^2 & v_1^2 & u_1v_1\\
1 & u_2 & v_2 & u_2^2 & v_2^2 & u_2v_2\\
\vdots & \vdots & \vdots & \vdots & \vdots & \vdots\\
1 & u_n & v_n & u_n^2 & v_n^2 & u_nv_n\\
\end{array}\right]$</center>
on fitte un polynôme de deux variables d'ordre deux.

On peut donc composer les régresseurs qu'on veut mettre dans la matrice $X$. A vous de voir ce qui correspond à votre problème : si j'ai bien compris vous avez bien deux variables pour prédire une sortie, vous êtes donc dans le cadre de l'exemple 2). Mais vous pouvez changer $X$ en rajoutant des ordres 3, ou 4, les interactions correspondantes etc.

Une fois la matrice $X$ définie, alors on peut estimer les coefficients $\beta$ par

<center>$\hat{\beta}=({}^tXX)^{-1}{}^tXY$</center>
ce qui se fait on ne peut plus simplement avec Matlab.
Enfin, on peut alors calculer l'estimation de $Y$ :

<center>$\hat{Y}=X\hat{\beta}$</center>

Est-ce suffisant ?
Dernière modification par Kuja le mardi 13 septembre 2005, 19:16, modifié 1 fois.
MB
Administrateur
Administrateur
Messages : 7729
Inscription : samedi 28 mai 2005, 14:23
Statut actuel : Enseignant

Message non lu par MB »

Kuja a écrit :Est-ce suffisant ?
Très bel effort et bonne utilisation du mode Latex. Merci Kuja.
MB. (rejoignez pCloud et bénéficiez de 10Go de stockage en ligne gratuits)
Pas d'aide en message privé. Merci de consulter ce sujet avant de poster votre premier message.
P.Fradin

Message non lu par P.Fradin »

Kuja a écrit :Une fois la matrice $X$ définie, alors on peut estimer les coefficients $\beta$ par

<center>$\hat{\beta}=(X'X)^{-1}X'Y$</center>
ce qui se fait on ne peut plus simplement avec Matlab.
Enfin, on peut alors calculer l'estimation de $Y$ :

<center>$\hat{Y}=X\hat{\beta}$</center>

Est-ce suffisant ?
Qu'est-ce que $X'$?
MB
Administrateur
Administrateur
Messages : 7729
Inscription : samedi 28 mai 2005, 14:23
Statut actuel : Enseignant

Message non lu par MB »

P.Fradin a écrit :Qu'est-ce que $X'$?
La transposée de $X$ non ?
MB. (rejoignez pCloud et bénéficiez de 10Go de stockage en ligne gratuits)
Pas d'aide en message privé. Merci de consulter ce sujet avant de poster votre premier message.
Kuja

Message non lu par Kuja »

Oui en effet, c'est bien la transposée de $X$.
Désolé d'avoir oublié de préciser la notation, c'est ma sale habitude d'utiliser Matlab (où la transposée s'écrit avec une apostrophe).
P.Fradin

Message non lu par P.Fradin »

Ne devrait-on pas la noter ${}^t\!X$?

Si c'est bien la transposée alors pour que la matrice ${}^t\!XX$ soit inversible, il y a une condition sur la matrice $X$, son rang doit être égal à son nombre de colonnes (ou noyau réduit à zéro).
Kuja

Message non lu par Kuja »

Bonsoir Pierre,

une fois de plus désolé pour la notation, la notation habituelle est bien sûr celle que vous proposez. Je vais donc éditer mon message de ce pas.
En ce qui concerne les conditions pour pouvoir inverser la matrice, elles sont en général vérifiées. En effet les régresseurs ne sont généralement pas choisis "au hasard", classiquement c'est une base (tronquée) de l'espace auquel appartient la fonction que l'on cherche à fitter.
P.Fradin

Message non lu par P.Fradin »

Kuja a écrit :Bonsoir Pierre,

En ce qui concerne les conditions pour pouvoir inverser la matrice, elles sont en général vérifiées. En effet les régresseurs ne sont généralement pas choisis "au hasard", classiquement c'est une base (tronquée) de l'espace auquel appartient la fonction que l'on cherche à fitter.
Merci du renseignement.

PS: moi c'est Patrick!
Kuja

Message non lu par Kuja »

Désolé pour le "Pierre", je ne sais pas d'où je l'ai sorti celui-là.
Bonne soirée Patrick :-)

Amicalement,
nirosis
Administrateur
Administrateur
Messages : 1803
Inscription : samedi 28 mai 2005, 14:48
Localisation : Orsay, France

Message non lu par nirosis »

Merci Kuja.

Cependant, mon problème est légèrelent différent, et je n'arrive pas à le modéliser correctement suivant le modèle classique.

J'ai $p$ fonctions du temps, $m_{1}(t),...,m_{p}(t)$ qui sont connues. Bien sûr, on a échantillonné le temps.

Le but est de montrer que ces $p$ fonctions sont corrélées et qu'en fait on peut les exprimer avec $q<p$ autres fonctions. C'est similaire à ce que pourrait donner une analyse en composantes principales, sauf que je veux qqchose de non-linéaire pour le coup. L'ACP permet d'écrire les fonctions $m_{i}(t)$ comme des combinaisons linéaires d'autres fonctions (en nombre le plus réduit possible). Le but est toujours de baisser le nombre de degré de liberté.

Je débute un peu sur ces problèmes de régression et d'analyse en composantes principales.
Je crois que la notion temporelle me trouble aussi, car les observations deviennent matrice et non plus vecteur.

ps : j'ai trouvé une thèse où des méthodes d'ACP non-linéaires sont utilisées. Je vais voir ce que ça donne. Lien
Kuja

Message non lu par Kuja »

Bonsoir nirosis,

effectivement étant données tes précisions ma réponse est complètement à l'ouest. Au moins j'aurais pu me faire plaisir sur ce forum :).
Pour en revenir à ton problème, je t'avouerais que j'ai un peu la flemme de me pencher sur l'ACP non-linéaire car je n'ai pas eu l'occasion de la rencontrer (je n'ai fait que de l'ACP classique). Donc si tu obtiens des résultats avec cette méthode, n'hésite pas à me (ou nous) les communiquer.
Cependant, il existe une méthode que j'ai déjà utilisée et qui marche plutôt bien, mais je ne sais pas si elle te suffira. C'est une méthode qui s'appelle ACE et qui sert à traiter les modèles linéaires généralisés (je doute quand même qu'on puisse faire beaucoup de choses pertinentes en sortant de ce cadre, car sinon les choses se compliquent très vite, c'est d'ailleurs un thème de recherche très actif en stat). En gros, étant donnés $p$ régresseurs $X_1,...,X_p$ et une sortie $Y$, la méthode ACE permet d'estimer les transformations $f_1,...,f_p$ et $g$ telles que

<center>$g(Y)=f_1(X_1)+...+f_p(X_p).$</center>
Si jamais tu ne trouves pas de références sur cette méthode (à condition qu'elle t'intéresse évidemment !), je devrais pouvoir remettre la main sur un article qui la présente.
nirosis
Administrateur
Administrateur
Messages : 1803
Inscription : samedi 28 mai 2005, 14:48
Localisation : Orsay, France

Message non lu par nirosis »

Salut, oui cela peut m'intéresser car ça généralise effectivement l'ACP classique.
En cherchant un peu j'ai trouvé cet article : Lien.
Je vais le lire un peu plus pour voir, mais je sais pas si c'est bien l'ACE dont tu parles.
Je veux bien que tu retrouves l'article qui détaille la méthode !

J'ai lu un peu les régressions non-linéaires : il y a des méthodes pour réduire à l'unidimensionnel... Et sinon ça utilise les réseaux de neurones, chose qui m'ennuie un peu (car je n'aime pas)... J'essairai peut etre un jour malgré tout.



ps : merci pour ton aide !
Biot

Détournement de topic, désolé :s

Message non lu par Biot »

Kuja a écrit :Bonjour,

en régression, le modélè qu'on considère est
<center>$Y=X\beta+\epsilon$</center>
où $Y={}^t(Y_1,...,Y_n)$ est le vecteur des $n$ observations, $X$ la matrice de taille $n\times p$ des $p$ régresseurs et $\epsilon$ le vecteur $n\times 1$ des erreurs. Dans la matrice des régresseurs $X$, on peut mettre ce qu'on veut. Quelques exemples :

1) Si on suppose que $Y$ est fonction d'une seule variable $x$ dont on a $n$ réalisations $(x_1,...,x_n)$ , alors en posant

<center>$X=\left[ \begin{array}{ccccc}
1 & x_1 & x_1^2 & ... & x_1^q\\
1 & x_2 & x_2^2 & ... & x_2^q\\
\vdots & \vdots & \vdots & \vdots & \vdots\\
1 & x_n & x_n^2 & ... & x_n^q\\
\end{array}\right]$</center>
on fitte un polynôme d'ordre $q$.

2) Si on suppose que $Y$ est fonction de deux variables $u$ et $v$ dont on a $n$ réalisations $(u_1,...,u_n)$ et $(v_1,...,v_n)$ , alors en posant

<center>$X=\left[ \begin{array}{cccccc}
1 & u_1 & v_1 & u_1^2 & v_1^2 & u_1v_1\\
1 & u_2 & v_2 & u_2^2 & v_2^2 & u_2v_2\\
\vdots & \vdots & \vdots & \vdots & \vdots & \vdots\\
1 & u_n & v_n & u_n^2 & v_n^2 & u_nv_n\\
\end{array}\right]$</center>
on fitte un polynôme de deux variables d'ordre deux.

On peut donc composer les régresseurs qu'on veut mettre dans la matrice $X$. A vous de voir ce qui correspond à votre problème : si j'ai bien compris vous avez bien deux variables pour prédire une sortie, vous êtes donc dans le cadre de l'exemple 2). Mais vous pouvez changer $X$ en rajoutant des ordres 3, ou 4, les interactions correspondantes etc.

Une fois la matrice $X$ définie, alors on peut estimer les coefficients $\beta$ par

<center>$\hat{\beta}=({}^tXX)^{-1}{}^tXY$</center>
ce qui se fait on ne peut plus simplement avec Matlab.
Enfin, on peut alors calculer l'estimation de $Y$ :

<center>$\hat{Y}=X\hat{\beta}$</center>

Est-ce suffisant ?
Bonjour, je me présente: Biot. Enfin c'est juste un pseudo pour ce forum quoi ^^

J'étais à la recherche d'informations pour calculer un polynome à deux inconnues qui me permette d'interpoler un nuage de points dans un espace à deux dimensions et mon ami google m'a finalement amené sur ce forum et plus particulièrement sur ce sujet.

Alors moi mon problème n'est pas celui posé par nirosis - et donc je m'excuse de détourner ainsi ce sujet - mais par contre la première réponse de Kuja s'applique très bien a ce que je recherche.

En fait j'appliquais déjà cette méthode avant de venir ici mais sans succès. Je pense qu'un bon point de départ pour trouver le polynome adéquat serait de savoir quel type de polynome on recherche, et c'est là que je coince. Comment savoir quel est le degré du polynome? Dois-je prendre le meme degré pour la variable x et pour la variable y? Quel produits intermédiaires ( xy, x²y, xy², x³y,...) dois-je inclure dedans?
J'ai essayé pas mal de combinaisons mais j'obtiens des formes qui n'ont rien a voir et meme parfois des amplitudes completement farfelues par rapport au points de départ. Le pire c'est que meme avec bcp de termes dans le polynome j'obtiens parfois des courbes étonnement plate alors que mon nuage de points est plutot chahuté.

Ma question est donc, étant donné un nuage de points, comment estimer quel type de polynome convient le mieux pour l'interpolation.


---

Autre chose, sur Matlab il existe une fonction interp2, qui marche à merveille pour l'exemple de matlab ainsi que pour un nuage de trois quatre points ... mais pas moyen de l'appliquer au problème qui m'interesse, matlab plante et/ou m'indique des erreurs (il faut dire que j'ai bcp bcp de points et que mon ordinateur n'est pas une bete de course).
Pensez vous qu'il soit intéressant d'insister avec cette fonction qui a l'air assez complexe par rapport à la méthode des moindres carrés que vous donnez et qui -elle- est toute simple à realiser dans matlab?

---

Encore autre chose, dans mes recherches (toujours avec google - je l'aime bien google, c'est mon pote) j'ai lu quelques articles sur la méthode du krigeage qui à été développée dans le cadre d'études géologique sur des gisements de minerais. Avez vous quelques connaissances là dessus? Cette solution fait appel à des notions de statistique et j'avoue que j'ai du mal a me replonger dedans.


voila, j'espère que je ne suis pas trop brouillon dans mes questions et que vous ne m'en voudrez pas de debarquer comme ça en plein sujet. Je ne verais pas d'objection à ce qu'on me dise de poster dans un nouveau sujet mais s'il vous plait prevenez moi avant de supprimer mon post. merci.


ps: aidez moi qd meme, j'ai pas envie de me faire battre par mon boss si je trouve pas ce maudit polynome :wink:
Biot

Message non lu par Biot »

Je me permets de faire un petit up au cas où ma question serait passée inaperçue :?
nirosis
Administrateur
Administrateur
Messages : 1803
Inscription : samedi 28 mai 2005, 14:48
Localisation : Orsay, France

Message non lu par nirosis »

Effectivement ton message est passé inaperçu pour moi...

Je crois qu'au départ il vaut mieux considérer toutes les puissances possibles. Ensuite s'il y a des coefficients très petits devant certains degrés, tu les enlèves et tu recommences avec un polynôme plus "adéquat" à ton nuage.

Le seuil pour éliminer un coefficient est déterminé par toi, tu peux prendre $1e-n$ avec un $n$ de ton choix... J'ai jamais testé reellement cette méthode, mais c'est ce qu'un prof m'avait dit.

Pour ta 2ème question, les fonctions matlab marchent bien encore faut-il comprendre ce qu'elles prennent et ce qu'elles renvoient... J'utilise souvent interp1 et c'est plutôt pas mal... Dis moi le message d'erreur pour savoir si c'est la mémoire qui sature ou autre chose... et combien as-tu de points ?
Biot

Message non lu par Biot »

Merci pour ta réponse.

Je me rends compte que mon problème serait peut etre plutot dû a une matrice mal continionnée lors du calcul de $\beta$, J'ai remarqué que matlab me renvoyait le warning suivant:

"Warning: Matrix is close to singular or badly scaled.
Results may be inaccurate. RCOND = 9.547710e-021."

Le problème doit venir lors du calcul de "$({}^tXX)^{-1}$" dans "$\hat{\beta}=({}^tXX)^{-1}{}^tXY$"


Pour ce qui est de la fonction interp2 de matlab, les problèmes que je rencontrais etaient principalement dûs au manque de mémoire de mon ordinateur ainsi qu'a la mise en forme des matrices input. Cela dit je préfère me concentrer sur la méthode des moindres carrés qui me permettrait d'obtenir un polynome et de ce fait trouver un maximum a ma fonction; la fonction interp2 ne donne finalement qu'une vue en 3D des données de départs, l'interpolation se limitant à relier les points par des plans.


sinon j'ai un echantillon d'a peu pres 1700 mesures; beaucoup ont des coordonnees x et y assez proches et des amplitudes assez cahotiques - je cumule les ennuis.
nirosis
Administrateur
Administrateur
Messages : 1803
Inscription : samedi 28 mai 2005, 14:48
Localisation : Orsay, France

Message non lu par nirosis »

Essaie de voir pourquoi cet inverse ne veut pas se faire... Affiche le determinant pour voir s'il est pas trop proche de 0...
Sinon calcul l'inverse d'une autre façon qu'avec la fonction matlab inv()...

Sinon tu peux essayer de remettre à l'échelle ton nuage de point en centrant tes variables (- moy / ecart type)... Ca peut aider à solutionner le pb de nuages "chaotiques".

1700 echantillons n'est pas énorme en soit cependant. Si un jour tu as l'occasion, passe ton script sur une machine plus puissante.