How can I regress Laguerre polynomials?
6 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
Marina
le 23 Fév 2014
Modifié(e) : John D'Errico
le 24 Fév 2014
Hi! I have to regress Laguerre polynomials, which have the following structure: L0=exp(-0.5*X); L1=exp(-0.5*X)*(1-X); L2=exp(-0.5*X)*(1-2*X+0.5*sqr(X)); LTot=aL0+bL1+cL2 where X is an MxN matrix Can you help me please?
0 commentaires
Réponse acceptée
John D'Errico
le 23 Fév 2014
Modifié(e) : John D'Errico
le 24 Fév 2014
I think you have a fundamental misunderstanding about this family of polynomials, and perhaps several misunderstandings.
(By the way, sqr is NOT a function in MATLAB, however, it tends to be easily misunderstood for sqrt, the square root. I don't know if you recognize that fact, and think that we should know that you mean x^2 when you write sqr(x), or if this is another mistake.)
Actually, what you have written are NOT Laguerre polynomials. The classic Laguerre weight function is exp(-x) here, such that they are orthogonal w.r.t. that weight function over the proper interval. Here the interval would be [0,inf). The polynomials are simply polynomials, thus 1, 1-x, etc. I'll discusss if you did even get the polynomials correct as we go.
For some reading, a good place to start is Abramowitz and Stegun, the true classic reference. A quick glance online might find other sources, such as wikipedia, or here .
Thus a Laguerre polynomial is a member of the polynomial family L_i(x), such that
int(L_i(x)*L_j(x)*exp(-x),[0,inf]) == d_(ij)
where d_(ij) is a delta function, equal to 1 when i==j, and equal to 0 when i~=j. The weight function is NOT part of the polynomial itself, although involved in the definition of the polynomial family.
So, are these Laguerre polynomials, for the given weight function? I recall having written a tool which generates the family of Laguerre polynomials. It resides in my sympoly toolbox.
orthpoly(0,'laguerre')
ans =
1
orthpoly(1,'laguerre')
ans =
1 - x
orthpoly(2,'laguerre')
ans =
1 - 2*x + 0.5*x^2
orthpoly(3,'laguerre')
ans =
1 - 3*x + 1.5*x^2 - 0.16667*x^3
These polynomials are such that they are orthogonal over the semi-infinite interval [0,inf), with respect to the weight function exp(-x). Lets test the orthogonality of the first few...
syms x
int(1*(1-x)*exp(-x),[0,inf])
ans =
0
int(1*(1)*exp(-x),[0,inf])
ans =
1
int((1-x)*(1-x)*exp(-x),[0,inf])
ans =
1
int((1-x)*(1-2*x+x^2/2)*exp(-x),[0,inf])
ans =
0
As you can see, these polynomials do appear to be orthogonal, at least for those I have tested. However, they are orthogonal for the weight function exp(-x). You have chosen a different family. If your weight function is exp(-.5*x) as you write it, then the polynomials are simply transformed from the classic ones I have written. Be careful though, as the family you have used is NOT orthogonal w.r.t. that weight function. The difference is perhaps subtle, but really quite important.
syms x
int(1*(1-x)*exp(-x/2),[0,inf])
ans =
-2
See that they were NOT orthogonal if we change the weight function. I'll give you the first few transformed Laguerre polynomials. The pattern should be obvious (I hope so, because I'm getting tired.) Thus, let
L0 = 1/sqrt(sym(2));
int(L0*L0*exp(-x/2),[0,inf])
ans =
1
L1 = (1 - x/2)/sqrt(sym(2));
int(L1*L1*exp(-x/2),[0,inf])
ans =
1
int(L0*L1*exp(-x/2),[0,inf])
ans =
0
I've chosen to normalize the polynomials by dividing by sqrt(2). This makes them truly an orthonormal set. You might have chosen as arbitrarily to make the constant coefficient always 1.
As far as regressing the polynomials, use backslash, creating the proper columns of your design matrix, as respectively
[ones(length(x),1)/sqrt(2),(1-x/2)/sqrt(2),(1-x+0.125*x^2)/sqrt(2),...]
Or use a tool from the stats toolbox for the regression. Be careful about the constant term though if you do, as regress presumes a unit constant term by default.
0 commentaires
Plus de réponses (0)
Voir également
Catégories
En savoir plus sur Polynomials dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!