This is something I’ve needed to do a lot recently on my current research project - calculate the area of polygons that lie within the boundaries of other polygons. For example, I wanted to find out the total area of arable land within each county/unitary authority in England and Wales. In a desktop GIS, this can be a fiddly, multi-step process and so it’s ripe for automation with some R code. The code snippet below shows how I went about this using functions from the tidyverse package and the new sf spatial data library. The post by “Matt SM” in this StackExchange thread was a useful starting point - http://gis.stackexchange.com/questions/140504/extracting-intersection-areas-in-r

#import required libraries
library(sf)
library(tidyverse)

#load county areas polygon shapefile 
counties <- st_read("counties_poly.shp") 

#load the arable land polygon data (extracted from the CEH LCM2007 1km land cover raster)
arable <- st_read("arable_poly.shp")

#before further processing, check the coordinate reference systems of both layers (must be same CRS)
#note: project to local coordinate system if current CRS in degrees
st_crs(counties)
st_crs(arable)

#run the intersect function, converting the output to a tibble in the process
int <- as_tibble(st_intersection(arable, counties))

#add in an area count column to the tibble (area of each arable poly in intersect layer)
int$areaArable <- st_area(int$geoms)

#plot the layers to visually check result of intersect
plot (arable$geometry, col='green')
plot(counties$geometry, add=T)
plot(int$geoms, col='red', add=T)

#group data by county area and calculate the total arable land area per county
#output as new tibble
tb_ArableByCounty <- int %>%
  group_by(County_UA) %>%
  summarise(areaArable = sum(areaArable))

#change data type of areaArable field to numeric (from 'S3: units' with m^2 suffix)
tb_ArableByCounty$areaArable <- as.numeric(tb_ArableByCounty$areaArable)

#join tibble to original county polygon shapefile and export as new shapefile
shp_out <- st_write(merge(counties, tb_ArableByCounty, by = 'County_UA'), "ArablebyCounty.shp")