badly written code may be slow (was hard to believe speed difference)
If you write R like C or Java, it will often be slow. R is a vector language, so learn to work with it rather than against it. For a truncated (standard) normal, try qnorm(runif(500000, a, b)) for suitable a and b (pnorm(lower) and pnorm(upper)?) And what are you going to do with 0.5m truncated normals that makes 70s significant? BTW, R has system.time() to time things.
On Tue, 6 Aug 2002, Jason Liao wrote:
First, I love R and am grateful to be using this free and extremely high quality software. Recently I have been comparing two algorithms and naturally I programmed in R first. It is so slow that I can almost feel its pain. So I decided to do a comparison with Java. To draw 500,0000 truncated normal, Java program takes 2 second and R takes 72 seconds. Not a computer science major, I find it hard to understand how R can be so slow.
It's not the language ....
R code:
normal.truncated <- function(mean, sd, lower, upper)
{
repeat
{
deviate = round( rnorm(1, mean, sd) );
if(deviate<=upper && deviate>=lower) break;
}
temp <- c(deviate+.5, deviate-.5, upper+.5, lower-.5);
result <- pnorm(temp, mean, sd);
prob <- (result[1] - result[2])/(result[3] - result[4]);
return(c(deviate, prob));
}
print(date());
n = 500000;
result = numeric(n);
for(i in 1:n) result[i] = normal.truncated(0,1, -10, 10)[1]
print(date());
Java Code:
package Mcdonald;
import VisualNumerics.math.*;
import java.util.*;
public final class normal_univariate_truncated
{
public double harvest, prob;
static Random ran_obj = new Random();
public normal_univariate_truncated(double mean, double var,
double lower, double upper)
{
double sd = Math.sqrt(var);
lower -= .5;
upper += .5;
double left = (lower - mean)/sd;
double right = (upper - mean)/sd;
double prob2 = Statistics.normalCdf(right) -
Statistics.normalCdf(left);
while(true)
{
harvest = mean + ran_obj.nextGaussian()*sd;
if(harvest>lower && harvest<upper) break;
}
harvest = Sfun.nearestInteger(harvest);
left = (harvest - .5 - mean)/sd;
right = (harvest + .5 - mean)/sd;
double prob1 = Statistics.normalCdf(right) -
Statistics.normalCdf(left);
prob = prob1/prob2;
}
public static void main(String[] useless)
{
System.out.println(new Date());
int n = 500000;
double[] result = new double[n];
for(int i=0; i<n; i++)
{
normal_univariate_truncated obj = new
normal_univariate_truncated(0, 1, -10, 10);
result[i] = obj.harvest;
}
System.out.println(new Date());
}
}
=====
Jason G. Liao, Ph.D.
Division of Biometrics
University of Medicine and Dentistry of New Jersey
335 George Street, Suite 2200
New Brunswick, NJ 08903-2688
phone (732) 235-8611, fax (732) 235-9777
http://www.geocities.com/jg_liao
__________________________________________________ -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
Brian D. Ripley, ripley at stats.ox.ac.uk Professor of Applied Statistics, http://www.stats.ox.ac.uk/~ripley/ University of Oxford, Tel: +44 1865 272861 (self) 1 South Parks Road, +44 1865 272860 (secr) Oxford OX1 3TG, UK Fax: +44 1865 272595 -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._