Village Occupation Modelling

Introduction

The goal of my project is to create a visualization of medieval village occupation based on data from test pit excavations.

Medieval rural settlement studies have traditionally focused on deserted settlements and overlooked currently occupied rural settlements. The deserted settlements could represent an atypical subset of medieval rural life (Taylor 1983), so currently occupied settlements are becoming an increasing research priority. Because archaeological data from currently occupied rural settlements is harder to access, test-pitting has been used to sample which portions of settlements were occupied over time (Jones and Page 2006, Gerrard and Aston 2007, Lewis 2014). Test-pitting locations are selected from residents who volunteer their gardens for two days of excavation with the goal of random distribution throughout the village (Lewis 2014).

Sampling of medieval settlements has been used to address a broad range of questions including continuity from the Roman period, the impact of the Black Death, how medieval villages formed, and whether medieval villagers had agency in the construction of their settlements. I would like to use modeling and visualization to create a representation of which parts of the village were occupied at a given time using test pit data.

This project will draw data from the University of Cambridge’s Access Cambridge Archaeology program of digging test pits in currently occupied medieval settlements. I selected this project’s data for this project primarily because of the strong commitment to publicly accessible data. Sherd count, total sherd weight, ware type, and context of recovery is all available on the publicly accessible pottery reports. Test pit location maps from each excavation season are also available. I selected the villages of Willingham and Peakirk

Load Packages

library(jsonlite)
library(maptools)
library(ggmap)
library(spatstat)
library(shp2graph)
library(plyr)
library(rgeos)
library(leaflet)
library(rgdal)
library(tmap)
library(ggplot2)

Geographic Setting

Test-pitting locations are generally located within eastern England. I decided to use Willingham and Peakirk for this project because more than 30 total test pit locations have been excavated for each village. The maps below illustrate the locations of these villages within the broader geography of England.

The shapefiles for the counties of England and Wales used to make these maps are available for download. Because accessing the data requires accepting terms and conditions, these maps load data from my computer rather than a public file source.

Map of Counties Highlighted within England

counties <- readShapeSpatial("~/Grad Year 1/Data Mapping/Final Project/countieswgs.shp")
counties[counties$G_NAME=='CAMBRIDGESHIRE AND ISLE OF ELY',] -> cambs
pt <- counties[counties$G_NAME == 'HUNTINGDONSHIRE AND PETERBOROUGH',]
proj4string(counties) <- CRS("+proj=longlat +datum=WGS84")
proj4string(cambs) <- CRS("+proj=longlat +datum=WGS84")
proj4string(pt) <- CRS("+proj=longlat +datum=WGS84")
tm_shape(counties) + tm_borders() + tm_shape(cambs) + tm_fill(col="dark green") +  tm_shape(pt) + tm_fill(col="dark green") + tm_grid(n.x=8, n.y = 8, col= "gray", labels.cex = 0.75, labels.col= "black", on.top=TRUE)

Map of the village locations and urban centers

Willingham is located 16km from Ely, while Peakirk is located about 6.5km from Peterborough. The map below shows the village locations with several urban centers shown purple in order to provide context.

Modern Willingham

Because medieval settlements comprised of a series of housing plots laid out along a street, the road layout will be an integral part to how I will calculate and visualize village occupation.

Modern housing developments with layouts such as the one near the intersection of Over road and Station road (highlighted in the inset) will interfere with how I plan to build my village occupation model.

library(ggmap)
get_map(location="Willingham, Cambridgeshire", source = "google", maptype="hybrid", zoom=15) -> map.will
get_map(location = c(lon = 0.05329, lat = 52.31223), color="color", source = "google", map = "roadmap", zoom=17) -> modern
ggmap(map.will) + inset_raster(modern, xmin=0.045, xmax=0.0525, ymin=52.306, ymax=52.310)

In order to eliminate modern roads from my village occupation models, I created a shapefile on QGIS by tracing over selected streets on the Open Layers Plugin. The shapefile only included the roads shown on the Ordnance Survey map of 1886. Nineteenth century Ordnance Survey maps are publicly accessible through archives such as this one. Although Willingham’s street plan changed and expanded between the end of the medieval period and 1886, my shapefile has excluded all roads that are obviously the result of 20th century housing developments.

Load Data

In addition to the roads shapefile mentioned above, I created two data files which will be used in the village occupation model visualization.

Test Pit Spatial Points

I created a shapefile with all excavated test pits shown on the map for excavations at Willingham in 2009 and 2013. I assigned numbers with the excavation year (9 or 13) followed by the number assigned to the test pit for the individual season.

tp <- readShapeSpatial("~/Grad Year 1/Data Mapping/Final Project/Village Models/WillinghamTPpointsall.shp")
tp <- tp[!is.na(tp$tpno),]

#save copy of data
w.tp <- tp

I have not been able to make the readOGR command work for GeoJSON, so the command line in my code involves reading the shapefile directly from my computer. I have the shapefile data available as GeoJSON here

Test Pit Occupation Data

I created a CSV compiling the data from the detailed breakdown of pottery recovered from each individual test pit. I will consider the question of how much pottery represents occupation (as opposed to sherds scattered in a field through manuring, for instance) later in this project. To build my model, I listed occupation based on whether the pottery report’s interpretion of each individual test pit represents occupation for each period. Some test pits produced more ambiguous data for inferring occupation, so I inputted my best estimate based on the wording of the report and the context data in the report. I listed the test pit as unoccupied if there was pottery present but pottery report described the location as a field. I included columns for both the total amount (sherd count and weight) of medieval pottery recovered and the amount (sherd count and weight) of pottery for each time period.

Because I will later be converting this data to an igraph object which will not accept empty rows, I decided to represent “yes” and “no” for overall occupation and occupation at each individual period with numeric values (1 for “yes” and -1 for “no”). Using numeric values now will permit me to enter 0 for all “NA” values in my data frame later using a single simple command. This is easier than entering values for empty spaces in character and numeric columns separately.

tp.df <- fromJSON("https://gist.githubusercontent.com/amw704/7e47c5cca57b63898f4a/raw/038f93769c3254462fc5154f53b00e2abd2437da/WillinghamTP.json")

#save copy of data
w.tp.df <- tp.df

The time periods were largely derived from the pottery report and are defined as follows: MS= Midddle Saxon: 650-850 AD; LS= Late Saxon: 850-1100 AD; M1= Medieval I: 1100-1200 AD; M2= Medieval II: 1200-1400 AD; M3= Medieval III: 1400-1550 AD

Roads Data

Using the process described in the “Modern Willingham” section, I created a shapefile of the late nineteenth century roads of Willingham:

roads <- readShapeSpatial("~/Grad Year 1/Data Mapping/Final Project/Village Models/WillinghamOSroads.shp")

#save copy of data
w.roads <- roads

I loaded the shapefile directly from my computer because I have not been able to load geojson with the readOGR function. The geojson file is available here

Map of Data Files

This leaflet map shows the roads shapefile traced over from the Ordnance Survey map and the test pit locations

leaflet() %>% addPolylines(data=roads, color="green") %>% addCircles(data=tp, color="blue", popup=tp$tpno) %>% addTiles()

Calculate Nearest Points Along Road to Test Pit Points

Medieval housing plots typically consisted of a toft and a croft and would have belonged to a household. The toft was the immediate area surrounding the medieval house and associated outbuildings (Browning and Higgins 2003). Uses of the toft included provisioning for animals, food storage for humans, quarrying for building materials, threshing floors for agricultural processing, and accommodation for elderly relatives (Astill 1988). If a test pit was interpreted as representing occupation (rather than use as a field), this interpretation implies that the area was part of a toft along the street. An occupied toft, in turn, implies that the portion of street in front of it was in use during medieval times.

Because village nucleation models generally attribute the development of the medieval street plan to the Late Saxon period at the earliest (Taylor 1983), I only used this modeling method for the periods post-dating the Late Saxon period. The test pits containing evidence of Middle Saxon occupation should not be expected to conform to the medieval or 19th century street plan because neither street plan would have existed during this time period. Only a single test pit produced Late Saxon evidence for Willingham, so the question of whether this method could be applied to the Late Saxon period is irrelevant for this village.

In order to make the calculation of the nearest point along the roads to each test pit, I converted the points to a ppp object and the roads to a psp object and used the project2segment function of the spatstat package.

#convert to ppp and psp objects
tp.ppp <- as(tp, "ppp")
roads.psp <- as(roads, "psp")

#calculation of nearest points along roads to test pits
tp.roads <- project2segment(tp.ppp, roads.psp)

tp.roads.df <- as.data.frame(tp.roads)

#save a copy of the data
w.tp.roads.df <- tp.roads.df

#create data frame with nearest points along roads and test pit number
tp.pts.only <- data.frame(x=tp.roads.df$Xproj.x, y=tp.roads.df$Xproj.y, tpno=tp.roads.df$Xproj.marks)

w.tp.pts.only <- w.tp.roads.df

#test plot to check that test pit points are now all along roads with original points in blue
plot(tp.pts.only$x, tp.pts.only$y)
plot(roads, col="red", add=TRUE)
plot(tp, col="blue", add=TRUE)

Convert Points and Roads to igraph object

I decided that the best method of connecting points along the roads while preserving the paths of the roads was to make an igraph object. With the “mapping.method” set to 2, the points2network function of the shp2graph function will add all of the test pit points as nodes to the network.

#the points2network function will only accept two columns of point data
pts.nodes <- tp.pts.only[,c("x", "y")]

#calculation of the length of the vector for the ea.prop value
roads.columns <- ncol(roads)

#conversion to network object
network <- points2network(roads, pts.nodes, mapping.method=2, ELComputed=TRUE, longlat=TRUE, ea.prop=(rep(c(0), length=roads.columns)))

#conversion to igraph object
graph <- nel2igraph(network[[1]], network[[2]], weight=network[[8]])

#test plot
plot(graph, vertex.size=3)

Extract Spatial Data to Join Nodes Representing Test Pits to the Test Pit Data Frame

Because the attribute “tpno” had to be removed to create the igraph object, it is necessary to extract the x and y values from all the nodes. The test pit numbers can then be joined to the data frame by matching the tp.pts.only data frame to the full list of nodes by x and y coordinates. When attributes are added to nodes in an igraph object, each value is matched in order to the corresponding node. For this reason, it is necessary to extract the list of nodes in order and create an attributes data frame based on the order of the nodes list.

#extract coordinates
V(graph)$x -> x
V(graph)$y -> y

#round by significant figures. I had trouble joining the nodes data to the test pit data because the precision of the two data frames differed
signif(x, digits=8) -> x
signif(y, digits=8) -> y

data.frame(x, y) -> nodes
signif(tp.pts.only$x, digits=8)-> tp.pts.only$x
signif(tp.pts.only$y, digits=8)-> tp.pts.only$y

#join test pit numbers to nodes list based on matching spatial data
join(nodes, tp.pts.only, match="first") -> nodes.tp.only

#join all test pit attributes from the json file to the ordered nodes data frame
join(nodes.tp.only, tp.df) -> nodes.df

#set all missing values to 0 so that R will perform the required calculations later
nodes.df[is.na(nodes.df)] <- 0

Add Attributes to Nodes Representing Test Pits

Now that the data has been joined based on the order of the points within the nodes list, it is possible to add attributes from the data frame to the nodes. For the purpose of building this model, I will only use the “occupation”" columns for each period based on the interpretation of each test pit within the pottery report. More detailed parameters of occupation based on the amount of pottery recovered would be possible by adding the pottery sherd and weight counts from each period to the nodes as attributes.

V(graph)$tpno <- nodes.df$tpno
V(graph)$occupation.ls <- nodes.df$occupation.ls
V(graph)$occupation.m1 <- nodes.df$occupation.m1
V(graph)$occupation.m2 <- nodes.df$occupation.m2
V(graph)$occupation.m3 <- nodes.df$occupation.m3

Get Shortest Paths

I could not identify a function on R that could connect all of the points representing occupied test pits and maintain the street pattern. It is important to preserve the street plan in the model since the medieval houses would have been laid out along the streets. Instead of connecting all of the points together, I used the get.all.shortest.paths function in igraph to calculate the shortest distance along the roads to the church (represented by number 999 in the data frame). The village church at Willingham has been in use since the Late Saxon period at the latest and would have been a central landmark in village life. It is therefore reasonable to assume that each occupied portion of the street (represented by the test pits) would have had a path to the church. The get.all.shortest.paths function highlights these paths for the church and then it is possible to extract these paths as a subgraph.

shortest <- get.all.shortest.paths(graph, from=V(graph)[tpno== "999"], to=V(graph)[occupation.m2=='1'])

#set shortest paths to red
for(p in shortest$res) {E(graph, path=p)$color <- "red"}

#testplot
plot(graph, vertex.size=3, vertex.label=NA)

#extract shortest paths by selecting red paths
short.graph <- subgraph.edges(graph=graph, eids=which(E(graph)$color=='red'), delete.vertices=TRUE)

Convert Subgraph to a Spatial Lines Data Frame

In order to use my outline of the village during the Medieval II period, I needed to convert the subgraph to a Spatial Lines Data Frame. I borrowed the code for this part from here

graph2_g<- short.graph
graph2 <-get.edgelist(graph2_g,names=F)
graph2<-data.frame(graph2)

graph2$Index <-c(1:nrow(graph2))
graph2$value <- E(graph2_g)$sum1
rownames(graph2)<-c(1:nrow(graph2))
colnames(graph2)<-c("X1","X2","Index")
graph2_shape <- apply(graph2,1,function(i) 
  Lines(Line(cbind(c(V(graph2_g)[i["X1"]]$x, V(graph2_g)[i["X2"]]$x)
                   ,c(V(graph2_g)[i["X1"]]$y, V(graph2_g)[i["X2"]]$y))),
        ID=as.character(i["Index"])))
graph2_shape <- SpatialLinesDataFrame(SpatialLines(graph2_shape), data=graph2,match.ID = F)

#test plot
plot(graph2_shape)

Visualization of Willingham 1200-1400

I could not directly plot the output of the previous step to map directly on leaflet (even after adding a projection), so I used writeOGR to create a shapefile from the Spatial Lines Data Frame. I then used readOGR to read the shapefile and plot it on leaflet.

This leaflet map represents the final output of the portion of Willingham occupied from 1200-1400 based on interpretations from the test pit data

willingham1200 <- graph2_shape
writeOGR(willingham1200, dsn=".", "willingham1200", driver="ESRI Shapefile", check_exists=TRUE, overwrite_layer=TRUE)
readShapeSpatial("willingham1200") -> willingham1200.shp
leaflet() %>% addPolylines(data=willingham1200.shp) %>% addTiles()

Repeat Steps for Medieval III period

I used the igraph object that I initially created and changed the get.all.shortest.paths calculation to select for Medieval III period (1400-1550AD) test pits instead. I kept all other aspects of the code the same, except for a subsitution of “willingham1400” for all instances of “willingham1200” in the names of the spatial lines data frame output, spatial lines length calculations, and shapefile output of the previous code. Because the code is nearly identical, I decided to hide most of it from my presentation.

shortest <- get.all.shortest.paths(graph, from=V(graph)[tpno== "999"], to=V(graph)[occupation.m3=='1'])
for(p in shortest$res) {E(graph, path=p)$color <- "blue"}
short.graph <- subgraph.edges(graph=graph, eids=which(E(graph)$color=='blue'), delete.vertices=TRUE)

Visualization of Willingham 1400-1550 AD

willingham1400 <- graph2_shape

writeOGR(willingham1400, dsn=".", "willingham1400", driver="ESRI Shapefile", check_exists=TRUE, overwrite_layer=TRUE)
readShapeSpatial("willingham1400") -> willingham1400.shp
leaflet() %>% addPolylines(data=willingham1400.shp) %>% addTiles()

Peakirk

Now that I managed to create visualizations of occupation estimates based on test pit data, I decided to repeat the process for another settlement, Peakirk. Peakirk is a significantly smaller settlement, which allowed for the Access Cambridge Archaeology team to place test pits closer together. Unlike Willingham, Peakirk contained significant evidence of occupation during the Late Saxon period, yet not a single test pit contained Middle Saxon evidence. The map below shows present-day Peakirk with the locations of testpit excavations.

roads <- readShapeSpatial("~/Grad Year 1/Data Mapping/Final Project/Village Models/Peakirk/Peakirkroads.shp")
tp <- readShapeSpatial("~/Grad Year 1/Data Mapping/Final Project/Village Models/Peakirk/PeakirkTPpoints.shp")
tp.df <- fromJSON("https://gist.githubusercontent.com/amw704/dfd8303be60ed114ab7b/raw/1ac983c1d684e64bb5ab446e5000941c3d9ff9e5/PeakirkTp.json")

#save copy of data
p.tp <- tp
p.tp.df <- tp.df
p.roads <- roads

As before, I could not use the readOGR function on my own computer, so I had to load the spatial data as shapefiles from my computer. The spatial data is available as GeoJSON. The roads data file is here and the test pit location data is here

This map contains both the test pit lcoations and the nineteenth century road plan.

leaflet() %>% addPolylines(data=roads, color="green") %>% addCircles(data=tp, color="blue", popup=tp$tpno) %>% addTiles()

#save a copy of data calculated from project2segment function
p.tp.roads.df <- tp.roads.df

Peakirk 1200-1400AD

peakirk1200 <- graph2_shape
proj4string(peakirk1200) <- CRS("+proj=longlat +datum=WGS84")
writeOGR(peakirk1200, dsn=".", "peakirk1200", driver="ESRI Shapefile", check_exists=TRUE, overwrite_layer=TRUE)
readShapeSpatial("peakirk1200") -> peakirk1200.shp
leaflet() %>% addPolylines(data=peakirk1200.shp) %>% addTiles()

Peakirk 1400-1550AD

peakirk1400 <- graph2_shape
proj4string(peakirk1400) <- CRS("+proj=longlat +datum=WGS84")
writeOGR(peakirk1400, dsn=".", "peakirk1400", driver="ESRI Shapefile", check_exists=TRUE, overwrite_layer=TRUE)
readShapeSpatial("peakirk1400") -> peakirk1400.shp
leaflet() %>% addPolylines(data=peakirk1400.shp) %>% addTiles()

Additional Criteria of Occupation

Based on published area excavations of medieval rural settlements, I decided to consider two additional criteria of whether test pit should be considered the result of occupation rather than manuring: distance from the street frontage and average sherd weight. By factoring these criteria into the determination of whether a test pit represents occupied area, I hope to develop a consistent method of addressing test pits for which the data is more ambiguous to interpret instead of using my best estimate. Adding these criteria should also make it easier for others to reproduce my results.

Relationship between distance from the street and pottery recovered

Fieldwalking surveys conducted during the Whittlewood research project demonstrated that the amount of pottery recovered generally decreased with the distance to the settlement (Dyer et al. 2005). I conducted my own spatial analysis of how much pottery per excavated square meter was recovered during the London Heathrow Terminal 5 excavations. The results of this spatial analysis also showed a decrease in the amount of pottery recovered with each successive interval of distance from the medieval settlement structures. A full explanation of these calculations is available here. In addition, the maps provided for the toft excavation at Antsey, Leicestershire show the boundary feature between the toft and croft at a distance of 40-60 meters from the road (Browning and Higgins 2003). The exact placement of toft boundaries could have been variable within a single settlement as at Raunds, NOrthamptonshire (Audouy and Chapman 2009), so it would be impossible to set a universal estimate of the distance between the outer toft boundary and the street frontage. On the other hand, this research suggests that test-pitting data will be a better indicator of occupation when placed within the plausible area of the toft.

I decided to see whether the test pit’s distance from the street frontage has an impact on the data from Willingham and Peakirk. I plotted the total weight of medieval pottery recovered on the y-axis with the distance from the street (calculated above using project2segment) on the x-axis. Test pits which produced evidence of Middle Saxon occupation were plotted in blue since this occupation happened before the medieval street frontages were laid out.

#use the saved output of project2segment to create a data frame of distance and test pit number
w.distance <- data.frame("tpno"=w.tp.roads.df$Xproj.marks, "distance" = w.tp.roads.df$d)

#combine with saved data frame of pottery recovered
join(w.distance, w.tp.df, by="tpno") -> w.distance2

#add column with village name so that these can be plotted with different markers
w.distance2$site <- "Willingham"

#repeat data frame preparation for Peakirk data
p.distance <- data.frame("tpno"=p.tp.roads.df$Xproj.marks, "distance" = p.tp.roads.df$d)
join(p.distance, p.tp.df, by="tpno") -> p.distance2
p.distance2$site <- "Peakirk"

#combine data from both settlements into a single data frame
rbind(w.distance2, p.distance2) -> total.dist
total.dist[is.na(total.dist)] <- 0
ggplot(total.dist[total.dist$weight > 0,], aes(x=distance, y=weight, color=occupation.ms, shape=site)) + geom_point()

This plot shows a general clustering of Peakirk test pits closer to the y-axis, which means that these tended to be excavated at locations close to the medieval street frontage. It also appears that Peakirk test pits produced more pottery than the Willingham test pits. This can be verified with a simple calculation:

#Willingham mean pottery weight per test pit, excluding test pits with 0 grams of pottery
#note: I excluded the two test pits representing Middle Saxon data since these areas would have been occupied before the development of the street frontage
mean(total.dist$weight[total.dist$site=="Willingham" & total.dist$weight > 0 & total.dist$occupation.ms == "-1"])

## [1] 33.63158

#Peakirk mean pottery weight per test pit, excluding test pits with 0 grams of pottery
#note: no Middle Saxon pottery was recovered at Peakirk, so it was not necessary to exclude Middle Saxon test pits
mean(total.dist$weight[total.dist$site=="Peakirk" & total.dist$weight > 0])

## [1] 58.11538

When pottery was recovered from a test pit, the average amount of pottery was greater in Peakirk than in Willingham. I decided to map how many test pits were within 50 meters of the street frontage for both villages. It is impossible to give an absolute distance of outer boundary of the tofts for an entire settlement, but fifty meters is a reasonable estimate based on the data cited above. If the distance of the test pits from the street frontage is impacting the amount of pottery recovered, I would expect to see a higher average of pottery recovered when the test pits further than 50 meters are excluded.

I first converted the projection of each data set to a projection that uses meters.

I used the rgeos package to calculate a 50 meter buffer around each set of roads. The sp::over function determined whether the test pits were within the buffer zone.

gBuffer(roads.wt, width=50) -> roads.wb
tp.wt$overlap <- over(tp.wt, roads.wb)
tp.wt$overlap[is.na(tp.wt$overlap)] <- 0
tp.wt[tp.wt$overlap=='1',] -> o.wt

gBuffer(roads.pt, width=50) -> roads.pb
tp.pt$overlap <- over(tp.pt, roads.pb)
tp.pt$overlap[is.na(tp.pt$overlap)] <- 0
tp.pt[tp.pt$overlap=='1',] -> o.pt

I created maps of the test pits and the 50 meter buffer zone. The test pits in red located within 50 meters of the street.

Willingham Test Pits within a 50 meter buffer

spTransform(o.wt, CRS("+proj=longlat +datum=WGS84")) -> o.wt2
spTransform(roads.wb, CRS("+proj=longlat +datum=WGS84")) -> roads.wb2

tm_shape(roads.wb2) + tm_borders(col="green") + tm_shape(w.tp) + tm_bubbles(size=0.25) + tm_shape(o.wt2) + tm_bubbles(size=0.25, col="red") + tm_grid(n.x=5, n.y = 5, col= "gray", labels.cex = 0.75, labels.col= "black", on.top=TRUE)

Peakirk Test Pits within a 50 meter buffer

spTransform(o.pt, CRS("+proj=longlat +datum=WGS84")) -> o.pt2
spTransform(roads.pb, CRS("+proj=longlat +datum=WGS84")) -> roads.pb2

tm_shape(roads.pb2) + tm_borders(col="green") + tm_shape(tp.pt) + tm_bubbles(size=0.25) + tm_shape(o.pt2) + tm_bubbles(size=0.25, col="red") + tm_grid(n.x=5, n.y = 5, col= "gray", labels.cex = 0.75, labels.col= "black", on.top=TRUE)

New average pottery recovery calculation

as.data.frame(o.wt2) -> o.col
join(o.col, w.tp.df) -> ow.df
ow.df$weight[is.na(ow.df$weight)] <- 0
mean(ow.df$weight[ow.df$weight > 0])

## [1] 37.625

The average amount of pottery recovered (for those test pits containing pottery) increased slightly when test pits outside the 50 meter buffer zone were excluded. It is worth noting the set of test pits which both had locations within the 50 meter buffer zone and contained pottery still includes test pits that were probably outside the area of occupation.

Only one Peakirk test pit (1234) was located outside of the 50 meter buffer zone and produced pottery. Because the total sherd weight was 37g, its exclusion should not significantly impact the average sherd weight for test pits within the 50 meter buffer zone. As a result, I decided not to recalculate the average pottery weight within the buffer zone for Peakirk.

Average Sherd Weight

In order to distinguish pottery deposition due to manuring from pottery deposition due to agriculture, the average weight of each sherd could be helpful factor. Plowing activity would have broken pottery sherds into smaller sherds over the course of the land’s use for agriculture (Lewis 2014). Paul Blinkhorn’s study on the mean sherd weight from the abandoned hamlet of West Cotton in Northamptonshire produced slightly lower average sherd weights in the yard areas (9.1-9.6 grams) when compared to the rest of the site (generally 10.1-10.5 grams). He attributes this difference to the higher rate of trampling within the yard areas in comparison to other parts of the site. He also notes that pottery discarded close to the site’s date of abandonment had a higher average sherd weight because these sherds would have been subject to less trampling (Blinkhorn 2010).

The average sherd weights for each individual buffer zone in the London Heathrow Terminal 5 spatial analysis produced a maximum average grams of pottery recovered per meter squared of 5.62. Rural occupation in this area continued until the mid-twentieth century, when it was obliterated for construction purposes. Perhaps the lower density of pottery recovery is due to the impact of trampling and waste disposal practices during the post-medieval period.

I created a column for the calculation of average sherd weight for both the Medieval II and Medieval III periods. I decided to use 4.0 grams as my tentative threshold of occupation as opposed to manuring.

w.tp.df$m2.sh.weight <- w.tp.df$m2.weight/w.tp.df$m2.sherd.co
w.tp.df$m3.sh.weight <- w.tp.df$m3.weight/w.tp.df$m3.sherd.co

p.tp.df$m2.sh.weight <- p.tp.df$m2.weight/p.tp.df$m2.sherd.co
p.tp.df$m3.sh.weight <- p.tp.df$m3.weight/p.tp.df$m3.sherd.co

Map of Willingham Test Pits and Average Sherd Weight

This map will show all Willingham test pits and highlight the test pits in red for which the average sherd weight is greater than 3.99 grams for the Medieval II or Medieval III period. The goal is simply to show how many test pits meet this requirement before factoring it into the village occupation model.

#add average sherd weight to shapefile and show average sherd weights greater than 3.99 in red. set NA values to 0
w.tp$m2.sh.weight <- w.tp.df$m2.sh.weight
w.tp$m3.sh.weight <- w.tp.df$m3.sh.weight
w.tp$m2.sh.weight[is.na(w.tp$m2.sh.weight)] <- 0
w.tp$m3.sh.weight[is.na(w.tp$m3.sh.weight)] <- 0

tm_shape(roads.wb2) + tm_borders(col="green") + tm_shape(w.tp) + tm_bubbles(size=0.25) + tm_shape(w.tp[w.tp$m2.sh.weight > 3.99 | w.tp$m3.sh.weight > 3.99,]) + tm_bubbles(size=0.25, col="red") + tm_grid(n.x=5, n.y = 5, col= "gray", labels.cex = 0.75, labels.col= "black", on.top=TRUE)

Map of Peakirk Test Pits and Average Sherd Weight

This map will show all test pits and highlight the test pits in red for which the average sherd weight is greater than 3.99 for the Medieval II or Medieval III period. The goal is simply to show how many test pits meet this requirement before factoring it into the village occupation model.

#add average sherd weight to shapefile and show average sherd weights greater than 3.99 in red. set NA values to 0
p.tp$m2.sh.weight <- p.tp.df$m2.sh.weight
p.tp$m3.sh.weight <- p.tp.df$m3.sh.weight
p.tp$m2.sh.weight[is.na(p.tp$m2.sh.weight)] <- 0
p.tp$m3.sh.weight[is.na(p.tp$m3.sh.weight)] <- 0

tm_shape(roads.pb2) + tm_borders(col="green") + tm_shape(p.tp) + tm_bubbles(size=0.25) + tm_shape(p.tp[p.tp$m2.sh.weight > 3.99 | p.tp$m3.sh.weight > 3.99,]) + tm_bubbles(size=0.25, col="red") + tm_grid(n.x=5, n.y = 5, col= "gray", labels.cex = 0.75, labels.col= "black", on.top=TRUE)

Create Another Village Occupation Model

I repeated the process of creating a village occupation model that factors in the 50 meter buffer zone and the sherd density. I added these as attributes to the data frame and only selected test pit data if they were located within a reasonable estimated of the toft area and had an average sherd weight of at least 3.99g

#join all test pit attributes from the json file to the ordered nodes data frame
w.tp.df$overlap <-tp.wt$overlap
join(nodes.tp.only, w.tp.df) -> nodes.df

#set all missing values to 0 so that R will perform the required calculations later
nodes.df[is.na(nodes.df)] <- 0

#add attributes to graph
V(graph)$tpno <- nodes.df$tpno
V(graph)$overlap <- nodes.df$overlap
V(graph)$m2.sh.weight <- nodes.df$m2.sh.weight
V(graph)$m3.sh.weight <- nodes.df$m3.sh.weight

I set the “get.all.shortest.paths” function to only use test pits within the 50 meter buffer zone and with a minimum average sherd weight of 4.0g

shortest.new <- get.all.shortest.paths(graph, from=V(graph)[tpno== "999"], to=V(graph)[overlap =="1" & m2.sh.weight > 3.99])

#set shortest paths to purple
for(p in shortest.new$res) {E(graph, path=p)$color <- "purple"}

#extract shortest paths by selecting purple paths
short.graph.new <- subgraph.edges(graph=graph, eids=which(E(graph)$color=='purple'), delete.vertices=TRUE)

I decided not to show the rest of the code from the rest of the script for building the model since it repeats steps from the previous example.

New Map of Willingham Occupation from 1200-1400AD

willingham1200.new <- graph2_shape
writeOGR(willingham1200.new, dsn=".", "willingham1200.new", driver="ESRI Shapefile", check_exists=TRUE, overwrite_layer=TRUE)
readShapeSpatial("willingham1200.new") -> willingham1200.new.shp
leaflet() %>% addPolylines(data=willingham1200.new.shp) %>% addTiles()

New Map of Willingham Occupation from 1400-1550 AD

I repeated the process, but I changed the code to select average sherd weights based on the Medieval III data instead of the Medieval II as before. Because the code is almost identical, I do not have it shown below.

willingham1400.new <- graph2_shape
writeOGR(willingham1400.new, dsn=".", "willingham1400.new", driver="ESRI Shapefile", check_exists=TRUE, overwrite_layer=TRUE)
readShapeSpatial("willingham1400.new") -> willingham1400.new.shp
leaflet() %>% addPolylines(data=willingham1400.new.shp) %>% addTiles()

The most obvious change in occupation between the two models is the absence of the extended north-south route in the 1400-1550 AD occupation model. One possible explanation for this is the existence of a mill owned by the landlord (the Bishop of Ely) near the area of the southern-most point of the 1200-1400 AD model, which is recorded as early as 1251. The mill was “under repair” in 1436. The absence of occupation evidence in the mill area during the 1400-1550 AD period could have been due to the disrepair of the mill, potentially for an extended period of time.

With these new criteria factored into the determination of whether a test pit represents occupation, the model selects significantly fewer test pits as occupied. The areas highlighted as occupied in each of the above maps are plausible areas of occupation based on the cultural resource management reports and the documentary evidence. The cultural resource management reports make it clear that more portions of the village were occupied than those represented by the village model. The village model does not include these because no test pits were dug in these portions of the village. Willingham is also currently one of the largest villages in Cambridgeshire, which probably makes sampling with test pits a more difficult task. It would be possible to add the data from the cultural resource management reports to the data frame in order to produce a more complete representation of village occupation based on all available excavated evidence in Willingham.

Conclusions

The development of this method of visualization of test pit data for a currently occupied rural settlements helped me think about the implications of sampling strategies. By linking together portions of the village interpreted as occupied during a given time, the model provides a different way of displaying village growth and decline as inferred by pottery data.

The biggest weakness of this modeling method is that the images of the Medieval II and Medieval III visualizations provide a static impression of occupation at Willingham over the course of two broad periods (200 years and 150 years respectively). This is drawback is more related to the often broad date ranges of the medieval pottery rather than the sampling methods. I also found that running the model on the Medieval I period produced an unrealistically large village occupation estimate. There are two probable reasons for this. First, occupation was inferred at two test pits dug in an area of Midddle Saxon settlement, where the test pits are located very far from the outlined roads. Additionally, the limited number of remaining test pits which produced Medieval I period evidence (8) were all located on the south end of the village at a significant distance to the church. If these test pits were occupied rather than part of an early field system, these could represent dispersed settlement as part of a gradual process of village nucleation.

It was useful to see the impact of adding the distance of the test pits to the streets and average sherd weight to the model. These criteria allow for a stronger argument that a specific test pit represents occupation, but the consideration of these factors also reduces the number of test pits that can contribute to creating the model. Because I have previously researched Willingham, I am already familiar with the results of area excavations conducted in advance of construction at various parts of the village. It would be interesting to take this modeling method further by integrating the cultural resource management data with the test pit data. This would create a data set with greater geographic coverage and a stronger basis for arguing when the various portions of the village were occupied.

Bibliography

Astill, Grenville. 1988. Rural Settlement: the Toft and the Croft. In The Countryside of Medieval England. Grenville Astill and Annie Grant, eds. Pp. 36-61. Oxford: Basil Blackwell Ltd.

Audouy, Michel and Andy Chapman. 2009. Raunds: The origin and growth of a midland village AD 450-1500. Excavations in Raunds, Northamptonshire 1977-87. Oxford: Oxbow Books.

Blinkhorn, Paul. 2010. The Saxon and medieval pottery. In West Cotton, Raunds: A study of medieval settlement dynamics AD 450-1450, Excavation of a deserted medieval hamlet in Northamptonshire 1985-89. Andy Chapman, ed. Pp. 259-334. Oxford: Oxbow Books.

Browning, Jennifer and Tim Higgins. 2003. Excavations of a Medieval Toft and Croft at Cropston Road, Anstey, Leicestershire. Transactions of the Leicestershire Archaeological and Historical Society 77: 65-81.

Dyer, Christopher, Richard Jones, Mark Page. 2005. Farming in Whittlewood AD43-1450: the Archaeological Perspective. In The Whittlewood Project: Medieval Settlements and Landscapes in the Whittlewood Area. Accessed online: archaeologydataservice.ac.uk.

Gerrard, Christopher M. and Michael Aston. 2007. The Shapwick Project, Somerset: A Rural Landscape Explored. Society for Medieval Archaeology. Leeds: Maney Publishing.

Jones, Richard and Mark Page. 2006. Medieval villages in an English landscape: Beginnings and ends. Oxford: Windgather Press.

Lewis, Carenza. 2014. The Power of Pits: Archaeology, Outreach and Research in Living Landscapes. In Living in the Landscape: Essays in Honour of Graeme Barker. Katherine Boyle, Ryan J. Rabett, and Chris O. Hunt, eds. McDonald Institute Monographs. Oxford: Oxbow Books.

Taylor, Christopher. 1983. Village and Farmstead: A History of Rural Settlement in England. London: George Philip.