[Rd] outer not applying a constant function

Martin Maechler maechler at stat.math.ethz.ch
Mon Mar 20 14:36:42 CET 2017


>>>>> Gebhardt, Albrecht <Albrecht.Gebhardt at aau.at>
>>>>>     on Sun, 19 Mar 2017 09:14:56 +0000 writes:

    > Hi,
    > the function outer can not apply a constant function as in the last line of the following example:

    >> xg <- 1:4
    >> yg <- 1:4
    >> fxyg <- outer(xg, yg, function(x,y) x*y)
    >> fconstg <- outer(xg, yg, function(x,y) 1.0)
    > Error in outer(xg, yg, function(x, y) 1) :
    > dims [product 16] do not match the length of object [1]

    > Of course there are simpler ways to construct a constant matrix, that is not my point.

    > It happens for me in the context of generating matrices of partial derivatives, and if on of these partial derivatives happens to be constant it fails.

    > So e.g this works:

    > library(Deriv)
    > f <- function(x,y) (x-1.5)*(y-1)*(x-1.8)+(y-1.9)^2*(x-1.1)^3
    > fx <- Deriv(f,"x")
    > fy <- Deriv(f,"y")
    > fxy <- Deriv(Deriv(f,"y"),"x")
    > fxx <- Deriv(Deriv(f,"x"),"x")
    > fyy <- Deriv(Deriv(f,"y"),"y")

    > fg   <- outer(xg,yg,f)
    > fxg  <- outer(xg,yg,fx)
    > fyg  <- outer(xg,yg,fy)
    > fxyg <- outer(xg,yg,fxy)
    > fxxg <- outer(xg,yg,fxx)
    > fyyg <- outer(xg,yg,fyy)

    > And with

    > f <- function(x,y) x+y

    > it stops working. Of course I can manually fix this for that special case, but thats not my point. I simply thought "outer" should be able to handle constant functions.

?outer   clearly states that  FUN  needs to be vectorized

but  function(x,y) 1    is not.

It is easy to solve by wrapping the function in Vectorize(.):

> x <- 1:3; y <- 1:4

> outer(x,y, function(x,y) 1)
Error in dim(robj) <- c(dX, dY) : 
  dims [product 12] do not match the length of object [1]

> outer(x,y, Vectorize(function(x,y) 1))
     [,1] [,2] [,3] [,4]
[1,]    1    1    1    1
[2,]    1    1    1    1
[3,]    1    1    1    1

----------------

So, your "should"  above must be read in the sense

  "It really would be convenient here and
   correspond to other "recycling" behavior of R"

and I agree with that, having experienced the same inconvenience
as you several times in the past.

outer() being a nice R-level function (i.e., no C speed up)
makes it easy to improve:

Adding something like the line

if(length(robj) == 1L) robj <- rep.int(robj, dX*dY)

before    dim(robj) <- c(dX, dY)   [which gave the error]

would solve the issue and not cost much (in the cases it is unneeded).

Or is this a bad idea?



More information about the R-devel mailing list