[SciPy-Dev] Improve Scipy's Stable PDF FFT Implementation

[Blair Azzopardi]

On Fri, 17 Aug 2018, 07:55 Andrea Karlova, <andrea.karlova at gmail.com> wrote:

> Blair,
>
> I think it may be good idea to finalise your note first and send it to
> academic journal,
> because it is their task to give you a feedback.
> Once you hear from them back, it is gonna be much easier to get approval
> for putting your preprint on ArXiv. Main idea behind ArXiv is to share
> preprints which are pretty much already in the publication process. For
> some journals, especially math ones, it takes years to get to the printed
> phase (my experience with good math journals is around 3 years in average).
>
> Regarding your note, I think you may need to add more text, literature
> overview, describe where your work is innovative, add section on future
> work and discussion and mainly formalise the body which you consider as a
> main results. Then just pick the journal and send it to them.
>
> It maybe also helpful to have particular code regarding the results on
> github, such that results from the paper can be easily replicable.
>
> Regarding the overall development: I have recently recieved email from the
> authors of the C-library libstable:
> https://www.researchgate.net/publication/317292451_libstable_Fast_Parallel_and_High-Precision_Computation_of_a-Stable_Distributions_in_R_CC_and_MATLAB
>
> They would be happy to collaborate on using their core library for SciPy,
> so that we would add wrappers, corrected math mistakes in the code with
> their help and they are happy to join and share their code with Scipy
> community. I think this is a robust longterm solution.
>
> Regards,
> Andrea Karlova
>
> On Tue, 14 Aug 2018 at 15:46, Blair Azzopardi <blairuk at gmail.com> wrote:
>
>> Hi
>>
>> I've written a short paper exploring optimising our Stable probability
>> density FFT implementation. In my original PR (7374) I introduced a method
>> by Mitnik. However, there is a better approach by Wang that significantly
>> improves accuracy. Wang's method can be generalized using standard
>> techniques improving accuracy even more but at the cost of CPU time. I felt
>> I should write to the dev mailing list as I'm not sure if implementing what
>> I feel is a generalized approach might not be established enough within the
>> academic community for inclusion in Scipy. Although, that said, Wang's
>> approach is cited by several other papers so shouldn't be a problem.
>>
>> The paper can be found in my github repository below. I'd like to submit
>> to Arxiv but can't without an endorsement; perhaps someone here might be
>> willing to help there.
>>
>>
>> https://github.com/bsdz/docs/blob/master/papers/Stable%20Probability%20Densities%20using%20FFTs%20and%20Newton-Cote.pdf
>>
>>
>> So, I plan to raise at least another PR to replace the FFT implementation
>> of Mitnik with Wang's. A couple of questions that come to mind:
>>
>> 1) Is it worth me just introducing the general solution that defaults to
>> Wang's?
>> 2) Is it worth placing the general FFT characteristic function code into
>> a separate module for approximating other integrals using FFT? If so, where
>> would be a good place?
>>
>> Blair
>>
>>
>> _______________________________________________
>> SciPy-Dev mailing list
>> SciPy-Dev at python.org
>> https://mail.python.org/mailman/listinfo/scipy-dev
>>
> _______________________________________________
> SciPy-Dev mailing list
> SciPy-Dev at python.org
> https://mail.python.org/mailman/listinfo/scipy-dev
>

Thanks for your comments. Just to clarify I didn't write the paper to be
innovative or for work towards a PhD. It's just a summary of known
available methods and a demonstration of how accurate they are compared to
what we have already and how we could improve on what we already have. I'm
happy to leave it in my GitHub repo.

The paper you linked to shows the same formulaes we already have although
they describe it as Nolan's method rather than my label 'Zolotarev'. I
added that implementation based on your concern that the existing DNI and
FFT methods weren't good enough if I recall correctly. Perhaps we should
rename "zolotarev" to "nolan" in the code.

One advantage to the approach described in that link you sent is it being
in raw C and utilising parallel processing. We could perhaps achieve
similar performance for our similar methods if we cythonize them; I'm not
sure if in Scipy we utilize mulithreading anywhere yet. Python has a few
issues with MT and GIL although doesn't stop using threads in compiled C
code. This might be an issue for portability.

The paper you link to doesn't really explore the FFT approach as much as
Wang has and if you read my investigation you'll see Wang's method achieves
much better accuracy than 10e-6; furthermore, by just cranking up the
cotes-newton degree we improve accuracy even more (10e-16) at the cost of
cpu time. It's well known that FFT method is better suited for alpha>1;
however, it's still useful for MLE of parameters.

Regarding the paper authors and libstable. In section 5 they state it's GPL
v3 and depends on GSL (gnu scientific library) so might require quite a lot
of effort for them to relax license to BSD compatible.

> <https://mail.python.org/mailman/listinfo/scipy-dev>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/scipy-dev/attachments/20180817/9113fb84/attachment-0001.html>