Estimating Volumes of Coastal Shell Midden Sites Using Geometric Solids

Abstract

Coastal shell midden deposits are a quintessential component of the archaeological record on the Pacific Northwest Coast. Despite their importance in informing the cultural and environmental histories of Indigenous peoples, research on shell middens has largely not sought to address the physical extent of these cultural deposits, which requires estimating shape, depth, and volume. Here, we present a new scalable geospatial model, designed to work with legacy survey data, for estimating midden volumes based on applying a regular geometric solid to sites with known extent and depth. We evaluate the accuracy of this technique using percussion core, total station, and lidar data from eight sites in Tseshaht territory on western Vancouver Island and three sites on the north coast of British Columbia (Canada). As part of the evaluation process of our results, we calculate uncertainty using subsurface core depth data and then compare generalized and modeled midden volume estimates. We demonstrate an accurate general model applied at the regional scale across a systematically surveyed landscape. This work presents the first landscape-scale measure of midden extents and volume within our study area, with relevance to historical ecology and settlement patterns.

Los concheros son un componente por excelencia del registro arqueológico en la Costa Noroeste del Pacífico. A pesar de su importancia para brindar información sobre las historias culturales y ambientales de los pueblos originarios, la investigación relativa a los concheros no ha buscado abordar la extensión física de estos depósitos culturales, motivo por el cual se requiere estimar su forma, profundidad y volumen. Aquí, presentamos un nuevo modelo geoespacial escalable, diseñado para trabajar con datos de estudios previos, a fin de estimar los volúmenes de estos depósitos, basados en la aplicación de una geometría regular del espacio a sitios con extensión y profundidad conocidas. Evaluamos la precisión de esta técnica utilizando datos de sondeos a percusión, estación total y lidar de ocho sitios en el territorio de Tseshaht en el oeste de la isla de Vancouver y tres sitios en la costa norte de la Columbia Británica (Canadá). Como parte del proceso de evaluación de nuestros resultados, calculamos la incertidumbre utilizando datos de profundidad del núcleo subsuperficial y luego comparamos estimaciones generalizadas y modeladas de los volúmenes de concheros. Demostramos un modelo general preciso aplicado a escala regional a través de un paisaje estudiado sistemáticamente. Este trabajo presenta la primera medida a escala de paisaje de las extensiones y volúmenes de concheros dentro de nuestra área de estudio, con relevancia para la ecología histórica y los patrones de asentamiento.

Shell middens are a globally distributed archaeological site type with considerable value to understanding the long-term history of Indigenous use of coastal ecosystems. Analysis of the constituents of shell midden deposits is a common research activity, but extrapolating data sampled from a small number of sites to the overall study area remains a methodological and interpretive challenge. Global research on shell middens is refining methods for linking evidence of settlement and resource use across a variety of settings (e.g., Astrup et al. 2019; Bailey and Hardy 2021; Litster et al. 2020; Parkington et al. 2020; Sanger et al. 2021). On the Northwest Coast of North America, shell middens are ubiquitous and have formative roles in cultural historical frameworks and zooarchaeological analyses (Ames and Maschner 1999; Cannon 2013). Although these cultural deposits are primarily composed of shell, bone, and other organic matter (Álvarez et al. 2011; Stein 1992; Waselkov 1987), they are also recognized as complex millennial-scale records of everyday life, including both the mundane and monumental with a variety of economic, social, political, and environmental significance (Álvarez et al. 2011; Blukis Onat 1985; Cannon 2000; Claassen 1991; Erlandson 2013; Grier 2014; Letham et al. 2020; Pluckhahn et al. 2015; Villagran 2014).

Coastal shell midden deposits can be large and topographically complex, which presents unique methodological challenges to survey and investigation (Claassen 1991; Lyman 1991; Stein 1992; Waselkov 1987). To date, little research has been undertaken that has explicitly looked at the spatial metrics and 3D shape of shell midden deposits. Most of the research that has been conducted has focused on recording the two-dimensional footprint of these features and sampling locations within them. Questions about the volume and shape of middens date back over 100 years (Nelson 1909; Spier 1916), and a few notable projects have attempted to estimate their physical extent and volume (Ceci 1984; Cook 1946; Letham et al. 2017; Randall 2015; Sorant and Shenkel 1984; Widmer 2014). Comparatively greater effort has been expended to explore the volume of terrestrial earthen mounds and other monumental features that are not shell bearing (Bernardini 2004; Blitz and Livingood 2004; Jeter 1984; Johnson 2015; Lacquement 2010; Magnani and Schroder 2015). Although the cultural contexts and formation process of earthen mounds differ substantially from shell middens, calculating the volume of these features involve similar principles and require engagement in the problem of modeling sediment shape and volume.

In this article, we present a new geometrically scaled model to estimate the volume and shape of coastal shell middens in sites on the Pacific Northwest Coast. This approach utilizes geographic information systems (GIS) and 3D modeling in combination with archaeological survey data to characterize and visualize site deposits. For validation, we use digital elevation models and subsurface testing data to create 3D models of shell middens, which can be used to calculate midden volumes for comparison to volumes estimated from a general shape. The implementation of these techniques allows for exploration of shell midden deposits at the landscape scale, and we present the first estimate of shell midden volumes along western Vancouver Island in British Columbia, Canada. This article concludes that a general model of midden shape can be created to estimate volume and spatial attributes for shell middens on the Northwest Coast and elsewhere. Additionally, we propose that this model is more accurate than previous efforts to employ geometric solids for modeling sediment volume.

STUDY AREA

The study area encompasses the Broken Group Islands, which are located in Barkley Sound, a marine embayment on western Vancouver Island. Covering approximately 100 km² and over 900 km of shoreline, the Broken Group Islands are a highly productive marine environment with more than 100 islands and islets fringed by rocky reefs, kelp forest habitat, and dense temperate rainforest (Figure 1). These islands are part of the traditional territory of the Tseshaht First Nation, a Nuu-chah-nulth group that occupied the region for millennia and now cooperatively manages this portion of Pacific Rim National Park Reserve with the federal Parks Canada agency (McMillan 1999; McMillan and St. Claire 2005).

Figure 1. The study area showing locations of shell midden sites. Lower inset map shows maximum depth of shell midden sites within the study area. The upper inset map shows the study area in relation to North America.

A considerable amount of archaeological survey and excavation has been conducted in the study area over the last 50 years (Haggarty and Inglis 1985; Hillis et al. 2020; McKechnie 2014; McMillan and St. Claire 1982, 2005; St. Claire 1975). In the early 1980s, the entirety of the Broken Group Islands was systematically surveyed as part of the process of establishing a cultural heritage baseline for Pacific Rim National Park Reserve (PRNPR; Haggarty and Inglis 1985). These pedestrian shoreline surveys in combination with additional cultural heritage monitoring and research over the subsequent 35 years have resulted in detailed site maps that record the locations and physical extent of shell-bearing midden deposits for 81 sites in the study area (Sumpter and St. Claire 2007, 2009; Sumpter and Suski-Armstrong 2005). In addition to this extensive history of archaeological work, this area also has a robust ethnographic and oral history record (Arima et al. 1991, 2000; Blenkinsop 1874; Drucker 1951; McKechnie 2015; Sapir and Swadesh 1939; St. Claire 1991). These data provide a strong foundation on which to consider research into the distribution of shell middens in the Broken Group Islands, and they are an extension of a multiyear collaborative research project between the Tseshaht First Nation, the University of Victoria, PRNPR, and the Bamfield Marine Sciences Centre (Tseshaht First Nation et al. 2017).

METHODS

We present a new process for modeling the volume of shell midden deposits using a geometric solid that is designed to be scaled across multiple sites where only midden location, extent, and maximum measured depth of cultural sediments are known. This process does not require extensive topographic surveying and can be applied at the regional level. Results are evaluated against primary data obtained from percussion coring at eight sites in the Broken Group Islands as well as three other sites on the northern British Columbia coast where similar surveys and subsurface testing have been conducted (Figure 2; Table 1). This validation process involves comparing volumes obtained from detailed archaeological investigation of these sites with the generalized geometric midden shape volumes.

Figure 2. The location of reference middens used to assess ellipsoidal volume estimates.

Table 1. Comparison of Volume Estimates.

^a Difference equals ellipsoid volume subtracted from reference volume.

^b Difference equals difference divided by reference volume (×100).

^c The absolute average error is 9% and is determined by summing the absolute value of each number in the %Diff column and dividing by the number of sites.

Midden Volume Estimation

To estimate shell midden volumes, we draw on the significant body of work employing the concept of a “geometric solid” (Shenkel 1986:201) as a representation of archaeological deposits based on idealized categories of geometric shapes (Ceci 1984; Cook 1946; Jeter 1984; Johnson 2015; Lacquement 2010). The geometric solid concept calculates volumes for geometric shapes and can be readily scaled across a range of sites without requiring highly specific topographic or photogrammetric data (Magnani and Schroder 2015). Provided that legacy mapping data for shell midden deposits are accurate, this approach can be used to estimate volumes for large numbers of sites with known boundaries. A geometric model will always be less accurate than a highly specific 3D model at an individual site (Álvarez et al. 2011), but as we describe below, our results improve on previous research by applying a new shape and developing a GIS-based method for linking the geometric solid to nonorthogonal footprints of shell-bearing archaeological deposits. The benefit of our approach is that it can also be applied without requiring lidar or other types of remote sensing data. Supplemental Script 1 includes the python code that was used to calculate shell midden volumes.

After a review of previous approaches, we observe that the generalized shape that seems to best represent shell-bearing midden deposits is a semitriaxial ellipsoid (Figure 3). This oblong and rounded shape is more specific and less symmetrical than other shapes previously used to estimate midden and mound volume such as cones (Sorant and Shenkel 1984), cylinders (Cook 1946), spheres (Treganza and Cook 1948), parabolas of revolution (Jeter 1984), pyramids (Shenkel 1986), and truncated “frustrum forms,” all of which do not correspond to the elongated crescent form of shell midden extents in our study area. The semiellipsoid shape can be scaled to the width, length, and depth of a midden with the assumption that its centroid corresponds to the deepest portion of the site and recedes in thickness to its horizontal margins at a curvilinear rate. We have observed that the deepest point in shell midden deposits tends to fall near the center, further indicating the suitability of this shape (Table 2). Given that coastal shell midden deposits tend to “smooth out” underlying landforms and shoreline topography, we truncate the ellipsoid into a half shape or “semitriaxial ellipsoid” to better represent a less topographically complex surface. A triaxial ellipsoid shape is defined by Equation 1:

(1)$$1 = \displaystyle{{x^2} \over {a^2}} + \displaystyle{{y^2} \over {b^2}} + \displaystyle{{z^2} \over {c^2}}$$

Figure 3. The semitriaxial ellipsoid: (a) the width of the ellipsoid as measured from the origin; (b) the length of the ellipsoid as measured from the origin; (c) the depth of ellipsoid as measured from the origin.

Table 2. Distance from Site Centroid to Deepest Core.

In Equation 1, x, y, and z represent the coordinates of the ellipsoid surface; a, b, and c are respectively the width, length, and depth of the ellipsoid as measured from the origin, which is placed at the site centroid (Figure 3). At each midden, the footprint is used to create a 1 × 1 m raster grid. Each cell in the raster is then converted to a point feature, thereby creating a regular grid of points within the midden, each of which represents 1 m² of space.

For each point feature in the grid, the x and y grid coordinates (UTM Zone 10 N) are recorded in the attribute table. As the UTM system records coordinates in meters north and east from an arbitrary origin, the raw coordinates enable volume calculations. The x and y distance from the midden centroid is used in Equation 2 (below) to calculate the height of the upper surface, or “skin,” of the ellipsoid at each point within the grid. Here, the variables represent the same attributes of the ellipsoid as Equation 1, with the addition of z, which is the vertical height (i.e., shell midden depth) of a point on the skin of the ellipsoid above a transverse plane drawn through the origin.

(2)$$z = \sqrt {\left({1-\displaystyle{{x^2} \over {a^2}}-\displaystyle{{y^2} \over {b^2}}} \right)\times c^2} $$

The a and b values in both equations are based on the width and length of a minimum bounding rectangle drawn around the site. The value of c is the maximum depth, as previously described. The height of the ellipsoid above the transverse plane can be calculated for any location on the ellipsoid's skin by solving for z. The UTM grid is oriented north, and the long axis of the shell midden deposit can have any orientation. For this reason, an additional adjustment must be made to the x and y values of each point to account for the difference between the orientation and the coordinate system (Figure 4). These x′ and y′ correction values are derived using Equations 3 and 4:

(3)$${x}^{\prime} = {-}x \times \cos ( \theta ) + y \times \sin ( \theta ) $$

(4)$${y}^{\prime} = {-}x \times \sin ( \theta ) -y \times \cos ( \theta ) $$

Figure 4. The steps of the volume estimation process: (a) shell midden site DfSh-31 with centroid, perimeter polyline, bounding rectangle, and raster grid; (b) the same site with points on a 1 m grid, showing angular adjustments to align the x′ and y′ calculations to the cardinal direction of the midden's long axis; θ represents the angular correction to the orientation of the midden; (c) midden points symbolized by depth with the ellipse layered over the midden; (d) core and inferred depth points placed on pseudo-transect lines at DfSh-31.

Each point created within the midden represents 1 m² of space, and it is assigned a height value representing the vertical distance from the transverse plane to the “skin” of the ellipsoid. The volume can be calculated by summing the z values of all the points on the grid within the ellipsoid. Recall that the volume contributed by each cell is z(m) ×  1 m ², or z(m ³). To calculate the volume of a full ellipsoid, this value would be doubled.

A benefit to an ellipsoidal shape is that it progressively thins toward the margins. However, in some cases, portions of a shell midden have irregular shapes that extend beyond the oval-like boundary of the ellipsoid. This is accounted for in our model through a spatial selection procedure that assigns the shallowest depth value found inside the ellipsoid footprint to those areas that extend beyond the ellipsoid perimeter. This manipulation is necessary to maintain the assumption that middens thin as they approach their edges. If left uncorrected, shell midden depth values begin to increase (rather than thin) as they extend beyond the ellipsoid boundary. This step effectively “clips” the footprint of the midden (Figure 4), like an oval cookie cutter passing through the shape of the midden where the area of the ellipsoid left inside the “cookie cutter” is what is used to calculate the midden volume. We believe this increases the realism and generality of the volume calculations when working with irregularly shaped middens while retaining the scalability and shape required for archaeological applications.

Our approach constrains the outline of the ellipsoid to the shell midden footprint and incorporates the desirable geometric qualities of both shapes. By preserving the unique shape of the midden footprint rather than using a standard geometric shape, we create a method that incorporates the individuality of discrete middens while also being applicable to a large number of sites—provided that the shell midden outlines are known. Previous methods that included this level of individuality (lidar or photogrammetry) cannot be as readily scaled, and they require unique calculations for each midden.

The shell midden settlements in this study include a wide variety of shapes, from relatively shallow sites with irregular boundaries to massive crescent-shaped mounded ridges hundreds of meters long and more than 5 m deep (McKechnie 2015). We did not encounter any sites in our study area that had footprints so irregular that they prevented the calculation of volume estimates using this method. That said, shell middens that do not fulfill the basic assumptions of being thickest near the highest point and gradually taper toward the edges would be less suitable for this method. Shell rings such as those found in the American Southeast would be an example for which this technique would not be appropriate (Thompson and Andrus 2011). However, for “stepped pyramid” mounds that do meet the stated assumptions, this technique would be appropriate (Pluckhahn et al. 2016). If a shell midden has been disturbed to the point where the footprint and depth cannot be measured (e.g., plowed-out middens), then this technique would not be applicable.

As can be discerned from the discussion above, we employed two foundational assumptions to develop our model. The first is that shell midden deposits are thickest nearest their centers and become thinner as they approach their margins, as variously observed through subsurface survey, mapping, and excavation on western Vancouver Island (Haggarty and Inglis 1985; McKechnie 2014). The second assumption is that shell midden deposits tend to occur on sloped shoreline landforms with a concave profile oriented toward the ocean. This is supported in the study area by exploratory data analysis of digitized midden footprints and lidar data, which shows that shell middens are predominately situated on sloped shoreline landforms. We employ the ellipsoid in our modeling because this shape captures both qualities of shell middens while retaining conservatism, in that a smoothed ellipsoid shape strikes a balance between overestimating volume on margins and underestimating volume taken up by irregular topography underlying shell midden deposits.

Location and Depth

We assembled digitized copies of original hand-drawn archaeological site maps from the Royal BC Museum survey records on file with the British Columbia Archaeology Branch (Haggarty and Inglis 1983, 1985). These data, including site reports and site locations, were used to roughly position site maps on the landscape. Site maps were loaded into a GIS and georectified using the scale bar on the original maps to resize the images to their correct size. Maps were then aligned to their real-world locations using digital elevation models and aerial imagery (Figure 5). Site maps for three shell midden sites out of the 81 in the Broken Group Island study area could not be located and were excluded from this analysis.

Figure 5. Site map georectification process: (a) the original site map for DfSh-31 (Haggarty and Inglis 1985); (b) the site map after being georectified; (c) the digitized midden footprint from the site map.

Once correctly positioned, midden extents were digitized by hand in ArcGIS Pro 2.7.1, creating a midden extent polygon layer. The majority of site maps for the study area were created using systematic surveying and mapping, and the internal accuracy of these records is very good, as described in Haggarty and Inglis (1985:55). Rubber-sheet georectification techniques were not employed because we found that these transformations introduced greater error than is present in the original maps. This is because most sites lack sufficient identifiable tie points to use in this process due to the heavy forest cover and the ephemeral nature of beach and coastline landscapes.

Midden depth data were compiled from a variety of sources including original survey reports, Parks Canada cultural resource management records, archaeological permit reports, and academic literature (Supplemental Table1). For each shell midden, the survey estimated or directly measured depth below surface of the basal cultural layer, which was recorded (Figure 1). Many sites in this area have had their depth measured. However, a smaller number have not been examined through subsurface depth testing, and for these locations, the reported depths are estimates based on observations of exposed stratigraphy at the site (Haggarty and Inglis 1985:55). For each site, the most specific depth measurement was used with preference given to depths generated from subsurface testing, and observational estimates were only used when a direct depth measurement was not available.

Eleven sites lack depth measurements and were estimated at a depth of 121 cm, which is the mean depth of sites within the study area. The use of mean depth is appropriate here because these sites include a wide variety of shapes and depths consistent with the larger sample population. We believe that the mean is the appropriate measure of central tendency to use here. For thoroughness, we also calculated volume estimates using median site depths and found that this change did not produce a significant change in volume (reducing the volume by about 1.5% on average).

Volume Estimate Assessment

To assess the accuracy of the volume estimates created from the semitriaxial ellipsoids, we compare these values to volume measurements created from excavation, percussion coring, augering, and total station mapping of midden sites. Although our sites have variable levels of sampling intensity, these comparisons allow us to assess how the model estimates perform as a measure of predicted shell midden volume. Sites used in this comparison were selected because they have many depth measurements that are well distributed across their extent.

Reference Volume Calculation

To assess the accuracy of the semitriaxial ellipsoid midden volume estimates, we developed a method of calculating shell midden volume from systematic subsurface surveys. A number of sites within the Broken Group Islands (n = 12) were cored (McKechnie 2014, 2015), and of these, eight sites were tested on a grid pattern. The data from this research provide detailed information about the basal depth of the shell middens and can be used to estimate volume. In order to create a more robust dataset to validate the model estimates, two middens that have been cored on a 10–20 m grid pattern in Prince Rupert Harbour (GbTo-34 and GbTo-70; Letham et al. 2017:7) and on the British Columbia Northern Coast (FhTj-1; McKechnie 2009) were added to our sample (Figure 2). From this point forward, we use the term “reference midden” to refer to these sites. A premise of this approach is that archaeological sites subjected to subsurface testing more accurately characterize sediment depths and volumes, which in turn can be compared with sites where only estimated depths and extents are known. We recognize that the total volume of a shell midden cannot be definitively known without complete excavation, and we consider this approach as coming closest to determining the total volume with existing data. The shell midden volumes from the Prince Rupert Harbour sites were used as reported, whereas volume calculations as described in the following paragraphs were performed on the raw data from the other nine sites. We have observed that a more thorough coring effort produces data that are better suited to measuring volume. As described by Martindale and colleagues (2009), percussion coring should obtain a large number of cores and be evenly distributed throughout the site to address variability in site depth, proportional to the size of the site and where smaller sites would not require as many depth measurements.

The first step in the calculation of reference midden volumes is to record the location of percussion cores and the maximum depth of cultural sediments from those cores in a GIS. Cores collected in the Broken Group Islands were analyzed at the Canadian Geologic Survey's Pacific Geoscience Centre, which produced high-resolution orthorectified photos. This imagery was used to determine the thickness of sands and gravels at the bottom of the cores, and this noncultural material was subtracted from total core length to derive the terminal depth of cultural sediments. The core depths and locations were then passed to the Inverse Distance Weighted (IDW) interpolation tool in ArcGIS Pro, creating a surface representing the bottom of the shell midden (Figure 6). IDW was determined to be the most appropriate method after comparing results using a variety of techniques (e.g., spline, kriging). This method has also been employed in other similar analyses (Johnson 2015; Monteleone 2007). For some sites, depth data from excavation units were also used to produce a more robust sample of midden depth measurements. High-resolution DEMs created from Total Station surveys (McKechnie 2014) and lidar data collected in 2018 allowed for detailed topographic surfaces to be created for each midden. Using these data, the volume of space between the top and bottom surfaces of the midden was determined using the ArcGIS Pro “Cut and Fill” tool (Figure 6; ESRI 2020). Although not as precise as conventional excavation, this volume calculation is derived from broadly spaced subsurface testing, and it offers an empirical basis for comparing with the ellipsoid estimation.

Figure 6. Visualizations of the reference and ellipsoid volumes for DfSh-31: (a) visualization of the top surface of DfSh-31; (b) visualization of the bottom surface; (c) the reference volume for this shell midden; (d) the semitriaxial ellipsoid volume of DfSh-31. The colored lines in B represent the location and depth of core tests used to derive the bottom surface.

A challenge in estimating volume from percussion core and topographic surface data is the need to use the inverse distance weighting (IDW) method to extrapolate the midden basal surface from known depths. As percussion cores surveys in the Broken Group Islands study area did not sample along the “edge” or outside site boundaries, we introduce additional “zero” data points following a “pseudo-transect approach” (Figure 4). Once these “pseudo-transect” lines were drawn, “assumed depth” points were placed at the site boundary and were given a depth of 0 cm. This is a variation of the technique that was used by Johnson (2015), who placed assumed zero-depth points at the vertices of midden boundaries. We modified Johnson's approach because the irregular midden shape can produce a high number of points and disproportionately influence IDW calculations, resulting in unrealistic thinning in large parts of the site. The pseudo-transect approach strikes a balance between edge points and core depths to be more reflective of the 3D shape and, ultimately, volume.

RESULTS

The integration of legacy survey data and geoarchaeological methods combined with geospatial analyses allows us to estimate site volumes in the Broken Group Islands. In addition, the digitization of sites in a GIS allowed us to determine that there is a high density of shell middens in the study area, with 3.93 sites per km² or roughly one site for every 430 m² of terrestrial terrain above low tide (Table 3). Following Mackie (2001), we additionally calculate the density per “coastal kilometer”—the Broken Group Islands encompass 299 km of linear coastline, which equates to 0.28 sites per km, or approximately one site for every 4 km of shoreline. Furthermore, a nearest-neighbor analysis shows that the average straight-line distance between sites is a mere 270 m.

Table 3. Shell midden Site Density Statistics Based on the Extent Shown in Figure 1.

Using the georectified site maps, the cumulative surface area of shell middens in the Broken Group Islands is estimated to be 107,253 m², which represents 0.5% of the total landmass (areas above low tide). The largest shell midden in the study has an area of 11,200 m², whereas the smallest encompasses 2.6 m², with the largest proportion tending to be smaller than 1,000 m². The distribution of depth, size, footprint area, and surface area are reported in Figure 7. The estimated cumulative volume of sites in the Broken Group Islands derived by the semitriaxial ellipsoid method is estimated to be 90,640 ± 2,719 m³.

Figure 7. Summary statistics for middens in the study area: (a) histograms of midden depth, (b) footprint area, (c) volume, and (d) surface area.

The uncertainty around midden volume calculations was determined by comparing the volumes generated from the ellipsoid models to the reference midden volumes across the nine sites (Table 1). This comparison suggests that, on average, the ellipsoid volumes have an absolute error of 9%. The net effect of under- and overestimation was 3% higher on average than the recorded volume. which we consider to be the margin of error for the volume calculations. Accordingly, the mean volume of a midden within the study area is 2,223 ± 67 m³, where the error margin is 3% of the volume. This increases to ± 200 m³ when the absolute error is considered. A large range in the size of shell middens with a skew toward large numbers of small shell middens results in a median value of 127 m³, which is considerably smaller than the mean volume. Such long-tailed distributions are not uncommon and were also reported in large-scale studies of shell middens north of the study area on western Vancouver Island and southern Haida Gwaii (Marshall 2006).

To assess the accuracy of reference midden volumes and how they vary with sampling intensity, we used percussion core data for the eight reference sites in the Broken Group Islands and a site on the North Coast (Table 4). We created a Python script in which 25% of cores were randomly removed, and then we recalculated volumes with the remaining depth measurements. To calculate coefficient of variation, the average standard deviation of the volumes from each site created by this process was divided by the average volume. This time-consuming process required over 24 hours to run and limited attempts at larger numbers of repetitions. However, the average variation is 13.6% of site volume, and site-specific variation can be seen in Table 4. Increasing the density of cores within a site increases the precision of the reference volume. We compared varying densities of cores per 100 m² of midden from our reference sites and found that the average coefficient of variation for a given density of cores increases as cores are removed in a linear fashion without any indication of what a suitable sample density might be (Supplemental Figure 1). This gives us a measure of confidence that we have sufficient core density in our samples for accurate calculation of volumes.

Table 4. Error for Reference Middens as Calculated by Adding and Removing Cores to Volume Calculations.

Note: Error is calculated as standard deviation of midden volume divided by the average of simulated site volumes.

The accuracy of the volume estimates was evaluated by comparing the volume as calculated for the semitriaxial ellipsoid against reference volumes calculated from systematic percussion coring at 11 sites (Table 1). This comparison shows that the semitriaxial ellipsoid–based estimates of midden volume overestimate reference midden volumes by 9% on average (absolute error), with divergence increasing between volumes predicted by the ellipsoid model and reference volumes at larger middens. For the three more intensively sampled reference middens outside the study area, the ellipsoid method tends to underestimate volume slightly (93%–94% of reference volume). For the Broken Group Islands, results are less consistent, likely reflecting the lower sampling resolution of percussion core depth measurement and the higher topographic complexity of the rocky intertidal shoreline in these wave-exposed island settings. At large sites (i.e., DfSh-4, DfSi-19, and DfSi-16), the ellipsoid estimate is between 11% and 22% larger than the reference volume, whereas for the other sites, the ellipsoid volume is less than 10% of the size of the reference volume. We observed that reference midden volume calculations perform better when more core data are available for a site because this increases the precision of interpolating midden basal surfaces. At the sites of DfSi-19 and DfSh-29, a relatively small number of cores were collected, whereas the other sites were more thoroughly surveyed, thereby producing a larger number of core tests. For instance, DfSh-79 is a “thin” but extensive defensive site located high on a steep-sided hill (McKechnie 2015), and it has a highly irregular footprint that does not conform to the usual crescent or oblong pattern observed in the study area. Its unique shape and elevated setting likely explain the discrepancy between the model and reference volumes. We believe that further testing on small sites and of midden footprint shapes that are more representative of those typically observed in the study area will refine our approach. Further attention to the topographic complexity of underlying landforms and bedrock exposures will also help refine estimates.

Overall, this method is considerably more accurate than previous applications of geometric solids for calculating archaeological midden and mound volume. To evaluate how our method compares to previous estimates, we applied previously proposed geometric solids to our reference middens and found that these methods vastly overestimated the sediment contained in middens in comparison to the reference volumes (Table 5). This overestimation can be anywhere from 177% to 926% on average, depending on the shape that is used. The shapes that were included in this analysis include cones (Sorant and Shenkel 1984), cylinders (Cook 1946), spheres (Treganza and Cook 1948), and pyramids (Shenkel 1986). This comparison to other geometric solids that have been used illustrates that semitriaxial ellipsoids provide substantially improved estimates of midden volumes over prior approaches.

Table 5. Comparison of Predicted Geometric Solid Models for Cacluating Site Volumes in Comparison with Measured Reference Volumes.

^a % Ref is calculated by dividing the predicted volume (m³) of a geometric solid in the column to the left by the corresponding value in the Reference Volume column and converting to a percentage.

DISCUSSION

The semitriaxial ellipsoid method used for estimating midden volume from geometric shapes improves the basis for evaluating the archaeological extent of shell middens on the Pacific Northwest Coast. The generic shape, scalable to site size, performs well in the absence of detailed topographic data or geophysical survey and is more appropriate for estimating sediment volumes than those used in previous efforts.

By employing a solid geometric shape, we recognize that we are sacrificing accuracy, as compared to systematic depth measurement across the extent of a midden, in favor of creating a process that can be applied to large numbers of sites—those with just two types of data: extent and depth estimates. Such information is routinely collected during archaeological survey efforts, which means that volume estimates can be generated from existing datasets. A benefit of applying a general solid is that it does not require high-resolution topographic data such as the kind produced by lidar or Total Mapping Station DEMs, which would be prohibitively time consuming and costly to acquire for all known sites in the densely forested Northwest Coast. For this same reason, geophysical surveys and photogrammetry are very difficult to conduct in this environment. In contrast, basic site maps have been created for thousands of recorded shell midden sites in British Columbia, and this method enables this extensive record to be analyzed with minimal additional fieldwork. The accuracy of the midden footprint and depth measurements will clearly influence the volume calculation. We do not recommend low-resolution data. Higher-resolution mapping produces more precise results, and we recommend using the best data available, but we suggest that the benefit rapidly diminishes beyond the submeter level. The resolutions of the input measurements are carried forward through the analysis to the final volume calculations. This method can be easily applied to very large collections of middens if they have been accurately mapped and they meet the model assumptions stated previously.

Through the application of GIS and Python scripting, volume calculations can be scaled across many sites simultaneously, with automated calculations allowing for examination of data on a landscape level. As discussed, empirical observations remain vital to refining this 3D shape characterization through systematic coring (Letham et al. 2017) and, perhaps, application of geophysical methods (Cariou et al. 2018; Kenady et al. 2018), although heavily forested sites on the Northwest Coast can be challenging (but see Urban and Carter [2021] for a recent successful application). Coring and auguring provide a good balance of time and effort (Cannon 2000; Duffield et al. 2020; Martindale et al. 2009; McKechnie 2014), and greater attention to this methodology will allow for refinement of accurate volume calculations in different cultural regions and coastal settings.

The ability to estimate midden volume accurately has significant potential for understanding shell middens and connecting results from excavated material to regional-scale processes. For instance, there is significant potential for shell midden volume estimates to be combined with radiocarbon dates and zooarchaeological data to better understand site formation processes, accumulation rates, and the density of animal remains at sites. This could shed light on the biomass represented in sites as well as the intensity of site use over time with implications for the modern use of marine resources (cf. Duffield et al. 2022; Reeder-Myers et al. 2022).

This is the first study that we are aware of to employ the semitriaxial ellipsoid as a model of shape to predict shell midden volume. The most similar shape that we have encountered is the “parabola of revolution” (Jeter 1984:92), which is superficially similar in shape but mathematically more similar to a cone. We argue that the ellipsoid better represents Northwest Coast shell midden deposits versus other shapes because the volume of an ellipsoid is derived from the measurement of three commonly described dimensions: length, width, and depth. Additionally, this shape is superior to a Cartesian model because the semirounded shape lacks linear planes, edges, and sharp corners. In contrast to the ellipsoid, the volume of spheres, cones, and similar shapes are determined using fewer axes of measurement. Additionally, our method constrains the ellipsoid to the footprint of the midden, which helps refine volumetric modeling of archaeological sediments. Together, these attributes of our model address several limitations of other approaches and enable more accurate volume estimates at scale.

Our analysis affirms that the Broken Group Islands were a densely inhabited coastal landscape and were substantially modified by the human activity represented in shell midden deposits. Over more than 5,000 years, Indigenous peoples living in this archipelago generated tens of thousands of cubic meters of sediment, including tens of millions of shellfish, fish, marine mammals, birds, and other organisms and cultural materials. The legacy of this long-term human activity is intimately connected to the historical ecology of this coastal region, but the extent and scale of this activity is not well recognized. We estimate the volume of archaeological shell midden sediments deposited in the Broken Group Islands (90,640 ± 2,719 m³) to be collectively equivalent to the volume estimates of famous monumental structures elsewhere in the Americas, such as the temples of Tikal (Webster and Kirker 1995) or the many mounds in the North American Southeast and Midwest (Shenkel 1986:214–215). These features in their totality represent the outcome of multiple generations of substantial labor and energy, especially considering the effort required to harvest materials in these deposits. Precontact modification to the landscape of the Broken Group Islands through the accumulation of sediments was not inconsequential, and we hope this study encourages archaeologists to consider ways that volume can enhance the interpretive and comparative potential of coastal archaeological landscapes.

CONCLUSION

Our geospatial analysis provides a novel method of estimating shell midden volume generated from readily available shape and depth information. Shell middens are often conspicuous deposits mapped with a considerable degree of accuracy and precision, but the combination of systematic subsurface testing, survey, and remotely sensed data to estimate site volume remains rare, despite significant potential on the Pacific Northwest Coast and elsewhere. We suggest the use and application of a semitriaxial ellipsoid model to estimate coastal shell midden volumes, and we compared these estimates with well-surveyed sites to establish the accuracy of this method. This method is scalable and can be applied to legacy data using an automated workflow. Estimating shell midden volumes provides interesting opportunities to better understand these complex anthropogenic features, which comprise major components of archaeological landscapes across the world. We hope that this approach helps spark discussions about characterizing the shape and extent of shell midden deposits and the historical ecology of human engagement in these different coastal regions.