Discussion :

[Natrio]

Bonjours à tous,

Mon instrument (en l'occurence un spectromètre) me délivre une valeur d'intensité en fonction de la longueur d'onde. Le spectre obtenu peut-être (doit-être!) approximé par une fonction de cette forme:

f(x) = e + f * cos(a + b*x + c*x*x + d*x*x*x)

Je cherche donc une méthode algorithmique me permettant de trouver automatiquement les coefficient a, b, c, d ,e et f. Leur ordre de grandeur peut-être fixé ou estimé a l'avance. Dans cette fonction, e et f correspondent respectivement à la valeur moyenne du signale et f son amplitude.

Quel méthode puis-je employer? J'ai regarder du coté des méthode de Newton ou quasi-newton mais je ne voie pas comment utiliser ces méthodes avec la fonction que je doit utiliser.

Je suis ouvert a toutes suggestions! Les méthode de newton me semblai bonne mais s'il y a une autre méthodes plus efficace ect...

Merci d'avance!

[pseudocode]

Regarde du coté des méthodes de régression non-lineaires (par exemple Levenberg-Marquardt)

[Natrio]

salut,

Suite tes conseils je me suis énormément servi de ce document :http://web.ics.purdue.edu/~kim497/ec...LM_handout.pdf pour mettre en place la méthode Levenberg-Marquardt sous Matlab.

Malheureusement cette méthode ne semble pas adapté a mon problème. Même en mettant des coefficients de départ très proche de ceux estimé manuellement la méthode trouve d'autre coefficients qui sont très mauvais.

Je cherche toujours et suis a l'affut de quelconque proposition!

[Natrio]

Je me demandais, es-ce qu'une régression au sens des moindre carré ne serait pas plus pertinente?

[pseudocode]

Je me demandais, es-ce qu'une régression au sens des moindre carré ne serait pas plus pertinente?

C'est justement ce que fait l'algorithme de Levenberg-Marquardt. C'est curieux que cela ne fonctionne pas dans ton cas.

[Natrio]

C'est curieux que cela ne fonctionne pas dans ton cas.

Bah vu ma maitrise des maths et de matlab il est très probable que mon algorithme soit mal codé...
Il me faut peut-être poster ça dans une autre section du forum (celle pour matlab?) Je poste ici au cas ou mais je déplace si nécessaire.

code:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
% initial guess for the parameters 
a0=0.00214892; b0=0.747274; c0=4.61907; d0=0.00667367; e0=0.0806972; f0=0.00109616;
a_est=a0; b_est=b0; c_est=c0; d_est=d0; e_est=e0; f_est=f0;
 
y_init = a0 + b0*cos(c0+(b0.*x)+(e0.*x.*x)+(f0.*x.*x.*x)); 
Ndata=length(y_input); 
Nparams=6;          % a b c d e and f are the parameters to be estimated 
n_iters=100;        % set # of iterations for the LM 
lamda=0.001;        % set an initial value of the damping factor for the LM
 
updateJ=1;
for it=1:n_iters 
    if updateJ==1   
        % Evaluate the Jacobian matrix at the current parameters (a_est, b_est) 
        J=zeros(Ndata,Nparams); 
        for i=1:length(x) 
            J(i,:)=[1, cos(c_est+d_est*x(i)+e_est*x(i)^2+f_est*x(i)^3), -b_est*sin(c_est+d_est*x(i)+e_est*x(i)^2+f_est*x(i)^3), -b_est*sin(c_est+d_est*x(i)+e_est*x(i)^2+f_est*x(i)^3)*x(i), -b_est*sin(c_est+d_est*x(i)+e_est*x(i)^2+f_est*x(i)^3)*x(i)^2, -b_est*sin(c_est+d_est*x(i)+e_est*x(i)^2+f_est*x(i)^3)*x(i)^3]; 
        end 
        % Evaluate the distance error at the current parameters     
        y_est = a_est + b_est*cos(c_est + (d_est.*x) + (e_est.*x.*x) + (f_est.*x.*x.*x)); 
        k=y_input-y_est; 
        % compute the approximated Hessian matrix, J’ is the transpose of J 
        H=J'*J;          
        if it==1              % the girst iteration : compute the total error 
            err=dot(k,k); 
        end 
    end 
 
    % Apply the damping factor to the Hessian matrix
    H_lm=H+(lamda*eye(Nparams,Nparams)); 
 
    % Compute the updated parameters 
    dp=-inv(H_lm)*(J'*k(:)); 
    a_lm=a_est+dp(1); 
    b_lm=b_est+dp(2); 
    c_lm=c_est+dp(3); 
    d_lm=d_est+dp(4); 
    e_lm=e_est+dp(5);
    f_lm=f_est+dp(6);  
 
    % Evaluate the total distance error at the updated parameters
 
    y_est_lm = a_lm + b_lm*cos(c_lm + (d_lm.*x) +(e_lm.*x.*x) + (f_lm.*x.*x.*x)); 
    k_lm=y_input-y_est_lm;     
    err_lm=dot(k_lm,k_lm); 
 
    % Ig the total distance error og the updated parameters is less than the previous one 
    % then makes the updated parameters to be the current parameters 
    % and decreases the value og the damping gactor   
 
    if err_lm<err 
        lamda=lamda/10; 
        a_est=a_lm; 
        b_est=b_lm;        
        c_est=c_lm;       
        d_est=d_lm;       
        e_est=e_lm;       
        f_est=f_lm;
        err=err_lm; 
        disp(err); 
        updateJ=1; 
    else % otherwise increases the value of the damping factor
        updateJ=0; 
        lamda=lamda*10; 
        a_est=a_lm;       
        b_est=b_lm;  
        c_est=c_lm; 
        d_est=d_lm;  
        e_est=e_lm;      
        f_est=f_lm;
        err=err_lm; 
        disp(err); 
    end  
end

Dison que les deux premier coefficient sont très faux (valeur moyenne et amplitude du signal) comparé aux autres appartenant au cosinus.

coefI =
e = 0.0021 f = 0.7473 a = 4.6191 b = 0.0067 c = 0.0807 d = 0.0011

coef =

e = 0.3203 f = 9.9780 a = 5.5640 b = 0.0846 c = 0.0774 d = 0.0007

coefI correspond au coefficient très proche de ceux déterminé manuellement, quand à coef, correspond au coefficient trouvé par la méthode Levenberg-Marquardt.

EDIT: Je vient de relire avec précision mon code, il s'avers que le calcule des dérivé partiel dans la matrice n'était pas tout a fait correcte. L'erreur rectifier permet une bien meilleurs estimation des coefficients.