Vectorised computation of sums of squares

Enrico Schumann

R, Matlab

## Review Status

In many applications the objective function is a sum of squares, that is for some vector $\epsilon$ of length $N$ we need to compute

(1)
\begin{align} \sum_i^{N} \epsilon_i^2\,. \end{align}

In vector-oriented languages like R and Matlab, matrix computations are often much more efficient than directly squaring the elements. Hence, replacing Equation (1) by

(2)
\begin{align} \epsilon'\epsilon\,, \end{align}

we get the same result (up to numerical precision), but generally the computation is much faster.

# R test

The following code was tested with R 2.8.1.

trials <- 10000
lenV <- 1000
aux <- rnorm(lenV)
#
t1 <- system.time(
for (i in 1:trials)
{
z1 <- sum(aux^2)
}
)
t2 <- system.time(
for (i in 1:trials)
{
z2 <- aux %*% aux
}
)
#
t3 <- system.time(
for (i in 1:trials)
{
z3 <- crossprod(aux)
}
)


Computing the inner product (method 2) already gives a considerable speedup, which rises to more than an order of magnitude with the function crossprod.

# Matlab test

The following code was tested with Matlab 2007a.

trials = 10000; lenV = 1000; aux = randn(lenV,1);

tic
for i = 1:trials
z1 = sum(aux .^ 2);
end
t1 = toc;

tic
for i = 1:trials
z2 = aux' * aux;
end
t2 = toc;

tic
for i = 1:trials
z3 = dot(aux,aux);
end
t3 = toc;


Here, the function dot (which is not a built-in function) is actually slowest.