Suggestions to speed up median() and has.na()

On Mon, 10 Apr 2006, Henrik Bengtsson wrote:

Hi,

I've got two suggestions how to speed up median() about 50%.  For all
iterative methods calling median() in the loops this has a major
impact.  The second suggestion will apply to other methods too.

I'm surprised this has a major impact -- in your example it takes much 
longer to generate the ten million numbers than to find the median.

Suggestion 1:
Replace the sort() calls with the .Internal(psort(x, partial)).   This
will avoid unnecessary overhead, especially an expensive second check
for NAs using any(is.na(x)).  Simple benchmarking with

x <- rnorm(10e6)
system.time(median(x))/system.time(median2(x))

where median2() is the function with the above replacements, gives
about 20-25% speed up.

There's something that seems a bit undesirable about having median() call 
the .Internal function for sort().

Suggestion 2:
Create a has.na(x) function to replace any(is.na(x)) that returns TRUE
as soon as a NA value is detected.  In the best case it returns after
the first index with TRUE, in the worst case it returns after the last
index N with FALSE.  The cost for is.na(x) is always O(N), and any()
in the best case O(1) and in the worst case O(N) (if any() is
implemented as I hope).  An has.na() function would be very useful
elsewhere too.

This sounds useful (though it has missed the deadline for 2.3.0).

It won't help if the typical case is no missing values, as you suggest, 
but it will be faster when there are missing values.

BTW, without having checked the source code, it looks like is.na() is
unnecessarily slow; is.na(sum(x)) is much faster than any(is.na(x)) on
a vector without NAs.  On the other hand, is.na(sum(x)) becomes
awfully slow if 'x' contains NAs.

I don't think  it is unnecessarily slow.  It has to dispatch methods and 
it has to make sure that matrix structure is preserved.  After that the 
code is just

     case REALSXP:
         for (i = 0; i < n; i++)
             LOGICAL(ans)[i] = ISNAN(REAL(x)[i]);
         break;

and it's hard to see how that can be improved. It does suggest that a 
faster anyNA() function would have to not be generic.