Using the R Standalone Math Library

Phil Spector

1  Introduction

One attractive feature of R is that the functions in its internal math library, including random number generators and probability functions, can be easily incorporated into C programs. At the SCF, the necessary library is available as described below. If you'd like to use this capability on your own computer, you need to download and configure R from source, and then change directories to $RHOME/src/nmath/standalone, and type make. This will produce a library called libRmath.so, which can be accessed through the command line in the usual way.

2  Using the Functions

While there is no formal documentation for using the R math functions, each of the functions accepts arguments directly corresponding to the arguments which the related R function accepts. Since there are no default parameters in C functions, it's important to provide all the arguments which the function needs. As a simple example, consider the quantiles of beta distribution. In R, the function that calculates these quantites is qbeta. First, the prototype for the function can be examined in the file $RHOME/include/Rmath.h:
% grep qbeta /usr/local/linux/R-2.11.1/include/Rmath.h
#define qbeta		Rf_qbeta
double	qbeta(double, double, double, int, int);

This shows that Rf_qbeta is the internal name which R uses to access this function - since we will be using the library in standalone mode, we don't need to worry about it. The function prototype for qbeta accepts two doubles and two ints. We can conveniently find what these arguments represent by looking at the source code for the qbeta function, along with the function's help page.
> qbeta
function (p, shape1, shape2, ncp = 0, lower.tail = TRUE, log.p = FALSE) 
{
    if (missing(ncp)) 
        .Internal(qbeta(p, shape1, shape2, lower.tail, log.p))
    else .Internal(qnbeta(p, shape1, shape2, ncp, lower.tail, 
        log.p))
}
<environment: namespace:stats>

This shows that, without a non-centrality parameter, we need to pass the desired probability, the two shape parameters, and two logical parameters. All logical parameters are internally passed as ints, with 0 meaning FALSE, and 1 meaning TRUE. The function prototype tells us that the function returns a double, so it's now clear how we need to call the function.

3  Calling qbeta from C

To insure that the necessary function prototypes are available to our program the Rmath.h file needs to be included. Since we are using the library in standalone mode, we need to define the variable MATHLIB_STANDALONE before including the header file. Once that's done, the rest of the program is quite straightforward. Suppose the following program is saved in the file beta.c:
#include <stdio.h>
#define MATHLIB_STANDALONE 
#include "Rmath.h"

main(){
   double shape1,shape2,prob;
   
   printf("Enter first shape parameter: ");
   scanf("%lf",&shape1);

   printf("Enter second shape parameter: ");
   scanf("%lf",&shape2);

   printf("Enter probability level: ");
   scanf("%lf",&prob);

   printf("Critical value is %lf\n",qbeta(prob,shape1,shape2,1,0));
}

To compile the program on the SCF Linux computers, use a compile line like the following:
cc -o beta beta.c -I/usr/local/linux/R-2.11.1/include -L/usr/local/linux/R-2.11.1/src/nmath/standalone/ -lRmath  -lm

A quick test shows that the program seems to work:
% ./beta
Enter first shape parameter: 1.5
Enter second shape parameter: 0.3
Enter probability level: .10
Critical value is 0.478398

This corresponds to the following call in R:
> qbeta(.10,1.5,.3)
[1] 0.4783981

4  Random Numbers

The following program illustrates calling an R random number generator from C. Note that the random number function (in this case rnorm) only returns one value, unlike the R function of the same name.
#include <stdio.h>
#include <time.h>
#define MATHLIB_STANDALONE
#include "Rmath.h"

main(){
   int i,n;
   double sum,x;
   time_t tt;
  
   tt = time(NULL);

/* the numbers passed to set_seed are arbitrary -- 
   use constants for a reproducible sequence       */

   set_seed(tt,77911);  

   printf("How many random normals to find the mean? ");
   scanf("%d",&n);

   sum = 0.;
   for(i=0;i<n;i++){
        sum += rnorm(0,1);      // mean, sd
        }
 
   printf("mean is %lf\n",sum/(double)n);
}




File translated from TEX by TTH, version 3.67.
On 8 Oct 2010, 11:06.