R/C Code

From QERM Wiki
Revision as of 21:00, 21 November 2012 by Chloe (Talk | contribs)
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
Jump to: navigation, search

Contents

Overview

Sometimes R users need to encounter C, for example, to maintain an R package, or to speed up existing functions. There are 2 options:

  • Write C code in a separate file, compile it, and call it from R
  • Use the package Rcpp to write C++ code within an R file

Note: this is one way to speed up your code, but don't forget about profiling your code (see R/Debugging) or writing it in a different language (e.g. Python or Java).

Writing C code separately

It is fairly straightforward to run C code from within R. Pdf10.png Here is how, with some easy-to-follow examples. More detailed information can be found in Writing R Extensions.

Tools needed

Mac/Linux

Everything should be already installed. You may need to install developer tools on Mac if they're not already there.

Windows

Install Rtools to get the necessary tools to compile code and build packages for R. Be careful if you already have Cygwin installed.

Simple example

From the command line:

R CMD SHLIB foo.c 

(you may need the whole path to R if it's not in your path)


To load what you just built to use in R:

dyn.load( file.path(path, paste("foo", .Platform$dynlib.ext, sep="")) )

Note: if you're on the mac and having issues with 32-bit vs. 64-bit (an error that says "mach-o, but wrong architecture"), you may need to specify the architecture: --arch=i386 for 32-bit and --arch=x86_64 for 64-bit, depending on which version of R.app you're running.

Using the Rcpp package

You can write C++ code within an R file and call it.

Simple example

Courtesy of Ingrid:

library(Rcpp)
library(inline)

#this is an example 

#here is my R function

mutate=function(Alleles,N,mu){  	
i=1
while(i<(2*N)){
prob=runif(1,0,1)    #First generate whether you are going to mutate
dir=runif(1,0,1)     #then decide which way to mutate (up or down) 50% chance either way 
if (prob<=mu&dir<=0.5){
Alleles[i]=Alleles[i]+1}
if (prob<=mu&dir>=0.5){
Alleles[i]=max(Alleles[i]-1,1)}
i=i+1
}
return(Alleles)
}

#here is the same function in Rcpp

mutateSrc = '
	using namespace Rcpp;
	
	RNGScope scope; 
	
	NumericVector Alleles(Al);
	int N = as<int>(Nind);
	double mu = as<double>(Mprob);
 	
	double prob;
	double dir;
	
	
	for (int i=0; i < 2*N; i++){
		prob=as<double>(runif(1,0,1));
		dir=as<double>(runif(1,0,1));		
		if(prob<mu&dir<0.5){	
			Alleles[i]=Alleles[i]+1;
		}else if (prob<mu&dir>0.5){
			Alleles[i]=Alleles[i]-1;
		}
				
	}

	return Alleles;
'

mutateRcpp = cxxfunction(signature(Al = "numeric",  Nind= "numeric", Mprob= "numeric"), body = mutateSrc, plugin = "Rcpp")

#here is a test with the same seed - same output and Rcpp is faster

set.seed(33455)

## test cases
system.time(
x <- mutate(rep(100,1000000),500000,0.9))

x[1:20]

set.seed(33455)

system.time(
y <- mutateRcpp(rep(100,1000000),500000,0.9))

y[1:20]
Personal tools