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.