R: calculating overlap between two standard deviational ellipses

Question

Aim

Calculating the area of overlap of two standard deviational ellipses in R and work out their intersection and union.

Tool to be used

calc_sde() from the 'aspace' package: calc_sde()

Question

Given some of the values returned by the function (see below), is it possible to calculate the area of the ellipses and the extent of their overlap? In particular, I would like to be able to work out the intersection and union between the two polygons. Not an expert in trigonometry and maths, so I am wondering which of the following values returned by the aboventioned function can be put to work and how:

CENTRE.x    
X-coordinate of the centre

CENTRE.y    
Y-coordinate of the centre

Sigma.x 
Half-length of axis along x-axis

Sigma.y 
Half-length of axis along y-axis


Theta   
Rotation angle in degrees

Eccentricity    
A measure of eccentricity (i.e., the flatness of the ellipse)

Area.sde    
Area of the SDE

TanTheta    
Trigonometric result

SinTheta    
Trigonometric result

CosTheta    
Trigonometric result

SinThetaCosTheta    
Trigonometric result

Sin2Theta   
Trigonometric result

Cos2Theta   
Trigonometric result

ThetaCorr   
Corrected theta angle for rotation of major axis from north

Answer 1

Lets make two sample objects, from the first and second five points in the sample data in the example for calc_sde:

z1 = calc_sde(id=1,calccentre=TRUE,weighted=FALSE, points=activities[1:5,], verbose=TRUE)
z2 = calc_sde(id=1,calccentre=TRUE,weighted=FALSE, points=activities[6:10,], verbose=TRUE)

That function is pretty appaling and I'd suggest you don't use it - it writes a file, it assigns objects in the global environment, and it prints that list out. All things that R functions really shouldn't do. Never mind. It returns a data frame that is the points of the ellipses, so we have the points in the ellipses now.

Using sf package functions we can create polygons:

library(sf)
s1 = st_polygon(list(as.matrix(z1[,2:3])))
s2 = st_polygon(list(as.matrix(z2[,2:3])))

And we can intersect them:

s12intersect = st_intersection(s1,s2)

plot(s1)
plot(s2, add=TRUE)
plot(s12intersect, add=TRUE, col="red")

And print the areas if we want that as well:

st_area(s1)
st_area(s2)
st_area(s12intersect)