Making an R package out of the C library.

[zoj613]

I'm wondering if there is enough interest in an equivalent R package.

The motivation would be that the existing R packages implementing the distribution do not also include functions for the CDF and density (as far as im aware). So maybe qpolyagamma and dpolyagamma would be useful for R users? Speed could also be a motivator?

The C library can be used to write a polyagamma sampler wrapper for any modern language that has support for C-extensions. The only thing it would need to implement separately is a bitgenerator to be used for random number generation and C implementations of random_gamma, random_standard_exponential and random_standard_normal. Both of the 2 requirements can be met fairly easily for an R extension by vendoring numpy's implementations found here and using the Xoroshiro128plus bitgenerator implementantion found in the C example. Nothing would need to modified in the C source of the library since the aforementioned functions are used via forward declaration and not included explicitly via header files.

If this R wrapper would be written, I think the following steps would be likely taken:

• create a seperate repo to house the R package.
• use this repo's C library (include and src folders) as a sub-module.
• add a distiributions.c file that vendors the random_gamma, random_standard_exponential and random_standard_normal implementations from numpy. Alternatively , other implementations could be used, but they would likely be slower.
• make sure the bitgen struct definition is available somewhere in the source:

     typedef struct bitgen {
       void *state;
       uint64_t (*next_uint64)(void *st);
       uint32_t (*next_uint32)(void *st);
       double (*next_double)(void *st);
       uint64_t (*next_raw)(void *st);
     } bitgen_t;

• Have a copy of bitgenerator implementation in the source, similar to what is in the C example.
• Write the r_wrapper.c code to make the C code available in R. See here for an example.
• Write the actual R implementations of ppolyagamma, dpolyagamma and qpolyagamma that call the C library functions. Due to R's design of random generation, we would need a polyagamma_set_seed utility function for setting the global seed.
• Optionally, this python interface which accepts a random_state parameter could be replicated fairly easily, but I dont think R users would like this kind of API.

NOTE: Something similar would be done if someone were to attempt to write a wrapper for another language like Julia.