[R-SIG-Finance] bessel functions

[Martin Maechler]

Hi Stefano,  et al.

I have been too busy to get this earlier,

>>>>> "Stefano" == stefano iacus <stefano.iacus at unimi.it>
>>>>>     on Thu, 14 Sep 2006 19:28:15 +0900 writes:

    Stefano> Hi, does anyone know about a reliable algorithm to
    Stefano> evaluate the exponentially rescaled bessel function
    Stefano> of first find, any order, real argument?

    Stefano> what I need is the equivalent of R function

    Stefano> besselI(x, nu, expon.scaled = TRUE)

    Stefano> for large real x (say 5000) and large nu (say 10 to
    Stefano> 30)

    Stefano> The R implementation does not satisfy my needs
    Stefano> because, for example,


    >> besselI(5539.947191,11,exp=TRUE)
    Stefano> [1] NaN

    Stefano> but if you ask Mathematica you get this number:
    Stefano> 0.00530180603357

    Stefano> To be more specific I need this to evaluate the
    Stefano> likelihood of the CIR model. I know there are many
    Stefano> different estimators for the CIR but I need to
    Stefano> evaluate the likelihood.

    Stefano> Googling around did not help me.

But  RSiteSearch("bessel")  does immediately help

The answer is *not* to use numerical recipes 
(they have ``source'' but have heavily restricted copyrights;
 and also very often, their code is far suboptimal) as 'mel'
indicated, but  for "special math functions" like these,
rather look at the GSL = GNU Scientific Library.

It has all kinds of Bessel functions and uses newer versions of
the algorithm than R, since
--- as we have the 'gsl' package on CRAN - thanks to Robin Hankin ---
you get the desired result from

 > library(gsl)
 > bessel_In_scaled(x= 5539.947191, n= 11)
 [1] 0.005301806


Martin Maechler