Created
November 19, 2012 15:26
-
-
Save hadley/4111256 to your computer and use it in GitHub Desktop.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
#include <Rcpp.h> | |
using namespace Rcpp; | |
// [[Rcpp::export]] | |
double vacc3a(double age, bool female, double ily){ | |
double p = 0.25 + 0.3 * 1 / (1 - exp(0.04 * age)) + 0.1 * ily; | |
p = p * (female ? 1.25 : 0.75); | |
p = std::max(p, 0.0); | |
p = std::min(p, 1.0); | |
return p; | |
} | |
// [[Rcpp::export]] | |
NumericVector vacc3(NumericVector age, LogicalVector female, NumericVector ily) { | |
int n = age.size(); | |
NumericVector out(n); | |
for(int i = 0; i < n; ++i) { | |
out[i] = vacc3a(age[i], female[i], ily[i]); | |
} | |
return out; | |
} | |
// [[Rcpp::export]] | |
NumericVector vacc4(NumericVector age, LogicalVector female, NumericVector ily) { | |
NumericVector p(age.size()); | |
p = 0.25 + 0.3 * 1 / (1 - exp(0.04 * age)) + 0.1 * ily; | |
p = p * ifelse(female, 1.25, 0.75); | |
p = pmax(0,p); | |
p = pmin(1,p); | |
return p; | |
} | |
// [[Rcpp::export]] | |
NumericVector vacc5(NumericVector age, LogicalVector female, NumericVector ily) { | |
NumericVector p(age.size()); | |
p = 0.25 + 0.3 * 1 / (1 - exp(0.04 * noNA(age))) + 0.1 * noNA(ily); | |
p = noNA(p) * ifelse(noNA(female), 1.25, 0.75); | |
p = pmax(0, noNA(p)); | |
p = pmin(1, noNA(p)); | |
return p; | |
} |
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
library(Rcpp) | |
library(microbenchmark) | |
vacc1a <- function(age, female, ily) { | |
p <- 0.25 + 0.3 * 1 / (1 - exp(0.04 * age)) + 0.1 * ily | |
p <- p * if (female) 1.25 else 0.75 | |
p <- max(0, p) | |
p <- min(1, p) | |
p | |
} | |
vacc1 <- function(age, female, ily) { | |
n <- length(age) | |
out <- numeric(n) | |
for (i in seq_len(n)) { | |
out[i] <- vacc1a(age[i], female[i], ily[i]) | |
} | |
out | |
} | |
vacc2 <- function(age, female, ily) { | |
p <- 0.25 + 0.3 * 1 / (1 - exp(0.04 * age)) + 0.1 * ily | |
p <- p * ifelse(female, 1.25, 0.75) | |
p <- pmax(0, p) | |
p <- pmin(1, p) | |
p | |
} | |
sourceCpp("vaccinate.cpp") | |
n <- 1000 | |
age <- rnorm(n, mean = 50, sd = 10) | |
female <- sample(c(T, F), n, rep = TRUE) | |
ily <- sample(c(T, F), n, prob = c(0.8, 0.2), rep = TRUE) | |
stopifnot( | |
all.equal(vacc1(age, female, ily), vacc2(age, female, ily)), | |
all.equal(vacc1(age, female, ily), vacc3(age, female, ily)), | |
all.equal(vacc1(age, female, ily), vacc4(age, female, ily)) | |
) | |
print(microbenchmark( | |
vacc1(age, female, ily), | |
vacc2(age, female, ily), | |
vacc3(age, female, ily), | |
vacc4(age, female, ily), | |
vacc5(age, female, ily) | |
)) |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment