'What is R's crossproduct function?

I feel stupid asking, but what is the intent of R's crossprod function with respect to vector inputs? I wanted to calculate the cross-product of two vectors in Euclidean space and mistakenly tried using crossprod .
One definition of the vector cross-product is N = |A|*|B|*sin(theta) where theta is the angle between the two vectors. (The direction of N is perpendicular to the A-B plane). Another way to calculate it is N = Ax*By - Ay*Bx .
base::crossprod clearly does not do this calculation, and in fact produces the vector dot-product of the two inputs sum(Ax*Bx, Ay*By).

So, I can easily write my own vectorxprod(A,B) function, but I can't figure out what crossprod is doing in general.

See also R - Compute Cross Product of Vectors (Physics)

r


Solution 1:[1]

Here is a short code snippet which works whenever the cross product makes sense: the 3D version returns a vector and the 2D version returns a scalar. If you just want simple code that gives the right answer without pulling in an external library, this is all you need.

# Compute the vector cross product between x and y, and return the components
# indexed by i.
CrossProduct3D <- function(x, y, i=1:3) {
  # Project inputs into 3D, since the cross product only makes sense in 3D.
  To3D <- function(x) head(c(x, rep(0, 3)), 3)
  x <- To3D(x)
  y <- To3D(y)

  # Indices should be treated cyclically (i.e., index 4 is "really" index 1, and
  # so on).  Index3D() lets us do that using R's convention of 1-based (rather
  # than 0-based) arrays.
  Index3D <- function(i) (i - 1) %% 3 + 1

  # The i'th component of the cross product is:
  # (x[i + 1] * y[i + 2]) - (x[i + 2] * y[i + 1])
  # as long as we treat the indices cyclically.
  return (x[Index3D(i + 1)] * y[Index3D(i + 2)] -
          x[Index3D(i + 2)] * y[Index3D(i + 1)])
}

CrossProduct2D <- function(x, y) CrossProduct3D(x, y, i=3)

Does it work?

Let's check a random example I found online:

> CrossProduct3D(c(3, -3, 1), c(4, 9, 2)) == c(-15, -2, 39)
[1] TRUE TRUE TRUE

Looks pretty good!

Why is this better than previous answers?

  • It's 3D (Carl's was 2D-only).
  • It's simple and idiomatic.
  • Nicely commented and formatted; hence, easy to understand

The downside is that the number '3' is hardcoded several times. Actually, this isn't such a bad thing, since it highlights the fact that the vector cross product is purely a 3D construct. Personally, I'd recommend ditching cross products entirely and learning Geometric Algebra instead. :)

Solution 2:[2]

The help ?crossprod explains it quite clearly. Take linear regression for example, for a model y = XB + e you want to find X'X, the product of X transpose and X. To get that, a simple call will suffice: crossprod(X) is the same as crossprod(X,X) is the same as t(X) %*% X. Also, crossprod can be used to find the dot product of two vectors.

Solution 3:[3]

In response to @Bryan Hanson's request, here's some Q&D code to calculate a vector crossproduct for two vectors in the plane. It's a bit messier to calculate the general 3-space vector crossproduct, or to extend to N-space. If you need those, you'll have to go to Wikipedia :-) .

crossvec <- function(x,y){
if(length(x)!=2 |length(y)!=2) stop('bad vectors')
 cv <-  x[1]*y[2]-x[2]*y[1]
return(invisible(cv))
}

Solution 4:[4]

Here is a minimalistic implementation for 3D vectors:

vector.cross <- function(a, b) {
    if(length(a)!=3 || length(b)!=3){
        stop("Cross product is only defined for 3D vectors.");
    }
    i1 <- c(2,3,1)
    i2 <- c(3,1,2)
    return (a[i1]*b[i2] - a[i2]*b[i1])
}

If you want to get the scalar "cross product" of 2D vectors u and v, you can do

vector.cross(c(u,0),c(v,0))[3]

Solution 5:[5]

There is a useful math operations package named pracma (https://rdrr.io/rforge/pracma/api/ or CRAN https://cran.r-project.org/web/packages/pracma/index.html). Easy to use and quick. The cross product is literally given by pracma::cross(x, y) for any two vectors.

Sources

This article follows the attribution requirements of Stack Overflow and is licensed under CC BY-SA 3.0.

Source: Stack Overflow

Solution Source
Solution 1
Solution 2 James Pringle
Solution 3 Carl Witthoft
Solution 4 Museful
Solution 5