Does anybody know some attempt to implement the Generalized Chi-squared distribution in Julia? (First I looked for it in Distributions.jl , but I have not found it there or anywhere else.)
Not that I’m aware of, but depending on what you want to use it for, a sampler should be easy to implement! If you want to go down the route of copying the R code to create a new Univariate Distribution I’m sure it would be much appreciated by the package maintainers
1 Like
Thanks. The R library seems to depend on C++ code that is not easy for me, and is GPL-licensed, so maybe it’s not a good choice. The Python library is MIT-licensed, but it seems also dependent on C++ and only implements the CDF, if I’m not mistaken. The Matlab library , which is also MIT-licensed, might be another option, I think.
1 Like
If the Matlab is MIT licensed, I’d definitely go for that! You can also try out the Matlab-to-Julia transpiler to make your job a little bit easier, although do note that it’s pretty rough and only gives a “first draft” that you’ll probably have to edit.
2 Likes
I have followed your suggestion:
JuliaStats:master ← heliosdrm:GChisq
opened 10:55PM - 27 Nov 23 UTC
This extends the collection of univariate continuous distributions with the Gene… ralized chi-squared. All required methods are added, and some of the recommended ones, as indicated in https://juliastats.org/Distributions.jl/stable/extends/#Univariate-Distribution
**Some notes:**
1. I started considering [Abhranil Das' code for Matlab](https://github.com/abhranildas/gx2) (MIT-licensed) as reference, although progressively departed from it, so that eventually there is little in the Julia code that can be attributed to him. I have kept a comment where I mention the inspiration to define the `cdf` function, using Davies' algorithm to integrate the characteristic function. But I don't know if a more explicit attribution in the code or in the license should be made.
2. There are several auxiliary functions that should not be part of the API, and I have encapsulated all of them in the module `GChisqComputations`, after the example of the Chernoff distribution.
3. That module includes two functions to calculate the characteristic function: `cf_explicit` that implements the explicit formula, and `cf_inherit` that composes the results of `cf(::Normal)` and `cf(::NoncentralChisq)`. Both are equivalent, and the code of the second is clearer, but `cf(::GeneralizedChisq)` calls the first one, because in a few tests that I have done, it is slightly (but not much) faster.
4. `cdf` and `pdf` are calculated by integration, using the Gil-Pelaez theorem (the algorithm for the CDF is given in Davies' paper, and the algorithm for the PDF has been indirectly derived from it). The integration is done using QuadGK.jl, which was already in the dependencies of Distributions.jl. The default tolerances of `QuadGK.quadgk` made the calculation of `cdf` at the median very slow, because there the integral is zero. Since in the ends of the distribution the integral is `±π/2`, which is in the order of magnitude of the unit, I have used a fixed absolute tolerance of `eps(one(T))` for the integration in the calculation of the CDF. I have not, however, changed the defaults for the calculation of the PDF.
5. `quantile` is calculated using `Distributions.quantile_newton`. However, there is no analytic solution to find the mode, which is the recommended starting point, so this is searched using the following strategy: (1) start at the mean of the distribution, and use it if it meets the convergence criteria; (2) otherwise define a bracket between that point and another at one standard deviation from the mean, in the direction towards the target probability (and if the bracket does not contain the target probability, extend it until it does); (3) bisect the bracket until one of the ends meets the convergence criteria, and then use that as the starting point. This seems to work fine, but it's a self-made algorithm, and perhaps there is a more efficient one that might be used instead.
2 Likes
