Program Accuracy, Limitations, and Versatility

The functions illustrated here are designed to have the look and feel of true mathematical functions. SPDF uses the Zolotarev transformation of the inverse Fourier transform of the characteristic function in the parameterization. The integrand is represented in the program symbolically so that transformation to the parameterization is exact. The SCDF and D1SPDF integrands are created by taking the integral and derivative of the SPDF integrand with respect to before proceeding with transform integration. The function, NIntegrate, gives a numerical approximation to the integral. It is important to remember that NIntegrate may not always come up with the correct answer. This is particularly true if the integrand has a sharp narrow spike. The integrand for the stable density can be problematic. The integrand is continuous and does not oscillate over the interval of integration; it begins with a value of 0, rises to a peak, , and falls to 0 at the other end of the interval. Unfortunately the peak for some values of , , and can be quite narrow and it is possible for NIntegrate to miss the peak with the default settings, resulting in an inaccurate approximation. Since the integrand is otherwise well behaved, it appears to be almost always possible to overcome this limitation and produce calculations to extreme accuracy. When the warning messages for NIntegrate are turned on, the first message will generally give a hint as to how to proceed. If the goal is high precision, the first strategy should be to increase the WorkingPrecision option. If machine precision is acceptable but the problem is that the peak is missed (most common with probability density function calculations), the best strategy is to increase the GaussPoints option. Even though the integrand does not contain a singularity, sometimes it will be necessary to increase the SingularityDepth option to catch the peak.

Here we could increase the recursion limit as suggested by the error message, but for these functions it is often a better idea to start with GaussPoints which sets a higher number of initial points, thus requiring less recursive subdivision. The difference between the calculations is small.

For extreme precision, start with working precision. As the output precision will be limited by the specification of the precision of the input values, it is simplest to use rational numbers. GaussPoints and/or MaxRecursion will also likely need to be increased. The default setting for PrecisionGoal in NIntegrate is WorkingPrecision minus 10 digits.

When the default NIntegrate parameters are used, the program may fail to give accurate results when is close to 1 or when is 1 and is close to 0. These limitations can often be overcome by changing the optional arguments to NIntegrate. With the warning messages for NIntegrate turned on, Mathematica will give notification that requested precision has not been achieved.

The package's Help Browser gives additional examples to help handle error messages generated by the other functions. It contains a section that demonstrates other methods of calculating stable distributions using Mathematica and includes the derivation of the Zolotarev integral. There is also an algorithm to set up an interpolation method for approximating the density function for maximum likelihood calculations.