[SciPy-User] orthonormal polynomials (cont.)

[josef.pktd at gmail.com]

On Tue, Jun 14, 2011 at 5:48 PM, Charles R Harris
<charlesr.harris at gmail.com> wrote:
>
>
> On Tue, Jun 14, 2011 at 3:15 PM, Charles R Harris
> <charlesr.harris at gmail.com> wrote:
>>
>>
>> On Tue, Jun 14, 2011 at 3:01 PM, <josef.pktd at gmail.com> wrote:
>>>
>>> On Tue, Jun 14, 2011 at 4:40 PM, Charles R Harris
>>> <charlesr.harris at gmail.com> wrote:
>>> >
>>> >
>>> > On Tue, Jun 14, 2011 at 12:10 PM, nicky van foreest
>>> > <vanforeest at gmail.com>
>>> > wrote:
>>> >>
>>> >> Hi,
>>> >>
>>> >> Without understanding the details... I recall from numerical recipes
>>> >> in C that Gram Schmidt is a very risky recipe. I don't know whether
>>> >> this advice also pertains to fitting polynomials, however,
>>>
>>> I read some warnings about the stability of Gram Schmidt, but the idea
>>> is that if I can choose the right weight function, then we should need
>>> only a few polynomials. So, I guess in that case numerical stability
>>> wouldn't be very relevant.
>>>
>>> >>
>>> >> Nicky
>>> >>
>>> >> On 14 June 2011 18:58,  <josef.pktd at gmail.com> wrote:
>>> >> > (I'm continuing the story with orthogonal polynomial density
>>> >> > estimation, and found a nice new paper
>>> >> > http://www.informaworld.com/smpp/content~db=all~content=a933669464 )
>>> >> >
>>> >> > Last time I managed to get orthonormal polynomials out of scipy with
>>> >> > weight 1, and it worked well for density estimation.
>>> >> >
>>> >> > Now, I would like to construct my own orthonormal polynomials for
>>> >> > arbitrary weights. (The weights represent a base density around
>>> >> > which
>>> >> > we make the polynomial expansion).
>>> >> >
>>> >> > The reference refers to Gram-Schmidt or Emmerson recurrence.
>>> >> >
>>> >> > Is there a reasonably easy way to get the polynomial coefficients
>>> >> > for
>>> >> > this with numscipython?
>>> >> >
>>> >
>>> > What do you mean by 'polynomial'? If you want the values of a set of
>>> > polynomials orthonormal on a given set of points, you want the 'q' in a
>>> > qr
>>> > factorization of a (row) weighted Vandermonde matrix.  However, I would
>>> > suggest using a weighted chebvander instead for numerical stability.
>>>
>>> Following your suggestion last time to use QR, I had figured out how
>>> to get the orthonormal basis for a given set of points.
>>> Now, I would like to get the functional version (not just for a given
>>> set of points), that is an orthonormal polynomial basis like Hermite,
>>> Legendre, Laguerre and Jacobi, only for any kind of weight function,
>>> where the weight function is chosen depending on the data.
>>>
>>
>> But in what basis? The columns of the inverse of 'R' in QR will give you
>> the orthonormal polynomials as series in whatever basis you used for the
>> columns of the pseudo-Vandermonde matrix.

with basis you mean hear for example the power series  (x**i,
i=0,1,..)  that's what they use, but there is also a reference to
using fourier polynomials which I haven't looked at for this case.

>>
>
> Example.
>
> In [1]: from numpy.polynomial.polynomial import polyvander
>
> In [2]: v = polyvander(linspace(-1, 1, 1000), 3)
>
> In [3]: around(inv(qr(v, mode='r'))*sqrt(1000./2), 5)
> Out[3]:
> array([[-0.70711,  0.     ,  0.79057, -0.     ],
>        [ 0.     ,  1.22352, -0.     , -2.80345],
>        [ 0.     ,  0.     , -2.36697,  0.     ],
>        [ 0.     ,  0.     ,  0.     ,  4.66309]])
>
> The columns are approx. the coefficients of the normalized Legendre
> functions as a power series up to a sign.

Looks interesting. It will take me a while to figure out what this
does, but I think I get the idea.

Thanks,

Josef

>
> Chuck
>
> _______________________________________________
> SciPy-User mailing list
> SciPy-User at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user
>
>