ENH: add Normal Inverse Gaussian (nig) distribution to scipy.stats (distributions.py).

The Normal Inverse Gaussian distribution gains quite a lot of momentum in the last years as a fat tail distribution (finance, Levy processes, extreme values, etc. see http://en.wikipedia.org/wiki/Normal-inverse_Gaussian_distribution)

I am mathematician and newbee in Python. Out of curiosity I coded with scipy a method to compute risk measures with nig distribution as base. Therefore I added the nig ppf to the library of continuous distribution functions in distributions.py. This is not difficult and I copied the code used for other distributions. However, being new to python, this may be not optimal.

I crossed checked with Matlab and R and so far it matches the result obtained there.

I am still thinking about a method to compute the cdf as well as the ppf since there are no close form for them and the performance are very poor in comparison with the equivalent in R or Matlab (This is central for our applications).

Being new to python, git and all the stuff I hope that I did the things at least halfway right.

~
ENH Normal Inverse Gaussian in scipy.stats

File distribution.py

According to the feedback I changed
1- removed spaces around ** (pep8 refuse that I remove \ for end of
lines)
2- removed m (or mu) as this is redundant with loc. Correct stats
accordingly
3- removed the redundant definition of logpdf and replace with
log(self._pdf(...))

Still open:
a) name (nig, norminv_gauss, norm_inv_gauss ?)
b) inclusion of init tests (I don't know which parameters I should
provide for for the distribution in the distribution list in
test_continuous_basic.py

General Problem
The computation of cdf and ppf is pretty slow in comparison to matlab or
R.
From my numerous tests, it seems that this has to do with the performance of the numerical integration
ENH: Normal Inverse Gaussian in scipy.stats

Following the discussion I changed:
1- change name to norminvgauss (seemed to be the consensus).
2- pdf is now 2 parameters
3- In comment added which parametrisation yields to the standard one.
4- In ./test/test_continuous_basic.py added norminvgauss to the list with parameters (2, -1)
ENH: modified norminvgauss and add genhyp to stats

Optimized the code for norminvgauss and add the generalized (proper)
Hyperbolic distribution.

Comments:
1- In the __call__ I initiate once for all some constants which are required to compute the pdf. This proves to be more efficient in terms of computation times. This is in particular the case for the genhyp which requires a computation of Bessel function for the constant (computation time almost divided by 2).

2- norminvgauss as well as genhyp does not provide very accurate values
for the cdf for high values due to complexity and singularity of the
    pdf. Therefore I also initiate for both in __call__ a variable x_opt
    which is the value of the mode of the distribution (maximum of the
    pdf). Since both distributions are unimodale, this mode theta is
    such that |theta - mean|<= sqrt(3) variance.
    Therefore I find the maximum of the logpdf on this interval using
    minimize_scalar. This is done only once, and afterwards the
    _cdf_single_call is computed via
    quad(_pdf, -infty,x) if x<=x_opt
    1-quad(_pdf, x,+infty) otherwise

    This provides significantly better computations of the cdf (and the
    ppf) in terms of accuracy as well as speed. The computation of the
    ppf in particular for high values does no longer fail and perform up
    to four times more quickly a result.

3- I rewrite the pdf using the exponential transform of the Bessel
function. This is better for the quad method and provides better results
as well as performance.

4- I crossed checked for several parameters with R and get the same
results. I will however pursue the test during the we.
MAINT: add norninvgauss to stats.distributions.__all__
MAINT: PEP8 and docstring fixes in stats.norminvgauss

Also:
- remove norminvgauss._logdf; it's the same as the generic implementation.
- create _compute_xopt function to remove code duplication.
- don't store ``fac1`` as attribute, it's cheap to compute.
MAINT: make calculating of x_opt in stats.norminvgauss work without state.