R help - Jun 2005 - A performance anomaly

I wrote a simple log likelihood (for the ordinary least squares (OLS)
model), in two ways. The first works out the likelihood. The second
merely calls the first, but after transforming the variance parameter,
so as to allow an unconstrained maximisation. So the second suffers a
slight cost for one exp() and then it pays the cost of calling the first.

I did performance measurement. One would expect the second version to
be slightly (very slightly) slower than the first. But I am finding
this is not the case! The second version is slightly faster. How can
this be?

On my machine (ibook @ 1.2 GHz, OS X "panther", R 2.1), I get:
> measurement(ols.lf1)
[1] 7.45
> measurement(ols.lf1.inlogs)
[1] 6.9

Here is the self-contained bug-demonstration code:

ols.lf1 <- function(theta, y, X) {
  beta <- theta[-1]
  sigma2 <- theta[1]
  if (sigma2 <= 0) return(NA)
  n <- nrow(X)
  e <- y - X%*%beta
  logl <- ((-n/2)*log(2*pi)) - ((n/2)*log(sigma2)) - ((t(e)%*%e)/(2*sigma2))
  return(-logl) # since optim() does minimisation by default.
}

# A variant: Shift to logs for sigma2 so I can do unconstrained maximisation --
ols.lf1.inlogs <- function(truetheta, y, X) {
  ols.lf1(c(exp(truetheta[1]), truetheta[-1]), y, X)
}

# We embark on one numerical experiment
set.seed(101)
X <- cbind(1, runif(20000))
theta.true <- c(2,4,6) # error variance = 2, intercept = 4, slope = 6.
y <- X %*% theta.true[-1] + sqrt(theta.true[1]) * rnorm(20000)

measurement <- function(f) {
  N <- 200
  cost <- system.time(
                      for (i in 1:N) {
                        j <- f(c(2,3,4), y, X)
                      }
                      )
  return(cost[1]*1000/N)                # cost per lf in milliseconds
}

measurement(ols.lf1)
measurement(ols.lf1.inlogs)

-- 
Ajay Shah                                                   Consultant
ajayshah at mayin.org                      Department of Economic Affairs
http://www.mayin.org/ajayshah           Ministry of Finance, New Delhi

Ajay Narottam Shah <ajayshah at mayin.org> writes:
> I wrote a simple log likelihood (for the ordinary least squares (OLS)
> model), in two ways. The first works out the likelihood. The second
> merely calls the first, but after transforming the variance parameter,
> so as to allow an unconstrained maximisation. So the second suffers a
> slight cost for one exp() and then it pays the cost of calling the first.
> 
> I did performance measurement. One would expect the second version to
> be slightly (very slightly) slower than the first. But I am finding
> this is not the case! The second version is slightly faster. How can
> this be?
> 
> On my machine (ibook @ 1.2 GHz, OS X "panther", R 2.1), I get:
> 
> > measurement(ols.lf1)
> [1] 7.45
> > measurement(ols.lf1.inlogs)
> [1] 6.9
> 
> Here is the self-contained bug-demonstration code:
What "bug"?

You're just not testing carefully enough. Timings are not perfectly
reproducible; you'll need multiple replications. I get
> replicate(50,c(measurement(ols.lf1.inlogs),measurement(ols.lf1)))
     [,1] [,2] [,3] [,4] [,5]  [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13]
[1,] 9.10 9.60  9.5 9.70 9.95 10.00  9.3 9.85 9.55  9.95 10.05 10.20  9.65
[2,] 9.75 9.65  9.2 9.15 9.60  9.55  9.8 9.85 9.75  9.70  9.35  9.65  9.60
     [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21] [,22] [,23] [,24] [,25]
[1,]  9.65  9.55   9.1  9.90 10.00  9.45  9.30   9.4   9.2 10.25  9.80 10.10
[2,]  9.95  9.90   9.7  9.55  9.55  9.80  9.95  10.6   9.3 10.30  9.55  9.65
     [,26] [,27] [,28] [,29] [,30] [,31] [,32] [,33] [,34] [,35] [,36] [,37]
[1,]   9.6  9.50  9.55  9.65  9.85  9.70  9.75   9.7  9.25  9.85   9.8  9.75
[2,]   9.7  9.75  9.10  9.75  9.40  9.55  9.95   9.9  9.55  9.55   9.3  9.60
     [,38] [,39] [,40] [,41] [,42] [,43] [,44] [,45] [,46] [,47] [,48] [,49]
[1,]   9.7  9.25  9.65  9.35  9.95  9.95  9.75  9.20  9.55  9.65  9.80  9.75
[2,]   9.9  9.65 10.20  9.20  9.65  9.50 10.25  9.75  9.90  9.10  9.25  9.05
     [,50]
[1,]  9.55
[2,]  9.40
> x <- .Last.value
> t.test(x[1,],x[2,])
        Welch Two Sample t-test

data:  x[1, ] and x[2, ]
t = 0.284, df = 96.202, p-value = 0.777
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -0.1018137  0.1358137
sample estimates:
mean of x mean of y
    9.663     9.646


-- 
   O__  ---- Peter Dalgaard             Blegdamsvej 3  
  c/ /'_ --- Dept. of Biostatistics     2200 Cph. N   
 (*) \(*) -- University of Copenhagen   Denmark      Ph: (+45) 35327918
~~~~~~~~~~ - (p.dalgaard at biostat.ku.dk)             FAX: (+45) 35327907