Vectorize!

Here is an email sent by one of my students a few days ago:

Do you know how to integrate a function with an  “if”?
For instance:
>X=rnorm(100)
>Femp=function(x){
+   return(sum(X<x))
+}
>integrate(Femp,0,1)$value

does not work.

My reply was that the fundamental reason it does not work is that integrate (or curve for instance) computes the function in several points of a grid by calling the function on a vector x and the comparison X<x does not make sense when both X and x are vectors. If one rewrites the code as

X=rnorm(100)
Femp=function(x){
   return((apply(as.matrix(outer(X,x,"-")<0),2,sum))}

the function can be integrated

> integrate(Femp,-2,2)
1.877735 with absolute error < 5.6e-05

However, this kind of syntactic gymnastic can be avoided by a simple call to Vectorize.

X=rnorm(100)
Femp=function(x){
   return(sum(X<x))}
> integrate(Vectorize(Femp),-2,2)
1.877735 with absolute error < 5.6e-05

2 Responses to “Vectorize!”

  1. I used a lot apply, but I did not know that Vectorize function… It looks great ! thanks….

    • Obviously, it is time-consuming:

      > system.time(for (t in 1:10^3) integrate(f,-1,1))
         user  system elapsed
        0.104   0.000   0.101
      > system.time(for (t in 1:10^3) integrate(Vectorize(f),-1,1))
         user  system elapsed
        0.416   0.004   0.418
      

      although a preliminary definition removes part of the time:

      > vf<-Vectorize(f,SIMPLIFY=TRUE) 
      > #does not work with FALSE
      > system.time(for (t in 1:10^3) integrate(vf,-1,1))
         user  system elapsed
        0.356   0.000   0.360
      

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s

This site uses Akismet to reduce spam. Learn how your comment data is processed.

%d bloggers like this: