Slow-moving, chronically destructive landslides are projected to grow in number globally in response to precipitation increases from climate change, and land disturbances from wildfire, mining and construction. In the Cincinnati and northern Kentucky metropolitan area, USA, landslides develop in colluvium that covers the steep slopes along the Ohio River and its tributaries. Here we quantify elevation changes in a slow-moving colluvial landslide over 14 years using county-wide lidar, uncrewed aerial vehicle (UAV) structure-from-motion (SfM) surveys and a UAV lidar survey. Because the technology and quality differ between surveys, the challenge was to calculate a threshold of detectable change for each survey combination. We introduce two methods; the first uses propagated elevation difference errors and the second back-calculates the individual survey errors. Thresholds of detection range from ±0.05 to ±0.20 m. Record rainfall in 2011 produced the largest vertical changes. Since then, the landslide toe has continued to deform, and the landslide has doubled its width by extending into a previously undisturbed slope. Although this study presents a technique to utilize older datasets in combination with modern surveys to monitor slow-moving landslides, it is broadly applicable to other studies where topographic data of differing quality are available.

Slow-moving landslides, which move at rates from just a few millimetres to several metres per year (Lacroix et al. 2020b) are chronically destructive and permanently damage property, infrastructure and agricultural land (Mansour et al. 2011; Nappo et al. 2019; Lacroix et al. 2020a). Although rarely deadly themselves, they can be precursors to fast-moving catastrophic landslides (Palmer 2017; Handwerger et al. 2019) and thus provide a valuable opportunity to study landslide processes prior to more rapid movement. The development of slow-moving landslides may result from perturbations to a steadily creeping slope (Chau 1999). Slow-moving landslides often remain inactive for years or decades, only to experience periods of rapid movement in response to precipitation, human disturbance or earthquakes (Lacroix et al. 2020b). With increases in precipitation owing to climate change, and land disturbance owing to wildfire, landslide activity is expected to increase in the USA (Leshchinsky et al. 2017; Mirus et al. 2017; Coe et al. 2018), and globally the number of landslides triggered by human activity such as construction and mining is increasing (Froude and Petley 2018). Within the USA alone, they are estimated to cost billions of dollars per year in economic losses (US Geological Survey 2005; Crawford 2014; Burns et al. 2017).

The purpose of the research we describe in this paper was to document and better understand decadal-scale spatial and temporal variations in the activity of a slow-moving landslide in the Cincinnati and northern Kentucky (USA) area, which necessarily included the development of an approach to integrate topographic data collected using different technologies over a period of years. Thus, this paper presents a method in a case study. Landslides throughout the area are well-known hazards but have not been studied intensively for nearly 30 years (Fleming and Taylor 1980; Haneberg 1991; Baum 1994; Fleming and Johnson 1994; Haneberg and Gokce 1994; Riestenberg 1994; Baum and Johnson 1996; Crawford 2012, 2014; Crawford and Bryson 2018; Glassmeyer and Shakoor 2021); perhaps because the landslides, although numerous and costly, are not deadly, and are treated largely as a continung maintenance problem. Thus, there are no published studies of the landslides incorporating modern remote sensing methods such as those we describe in this paper. The previous studies were also typically limited to a year or two in duration, so they did not provide information about decadal patterns of activity. This paper provides the first documentation of decadal Cincinnati and northern Kentucky area landslide movement based upon modern remote sensing techniques.

We describe a method to quantify the threshold of detectable change between digital elevation models (DEMs) produced using various remote sensing technologies. Here we use the US Geological Survey (USGS) definition of a DEM, which is a bare-earth surface that excludes vegetation, buildings and other surface objects. This is equivalent to a digital terrain model (DTM), as opposed to a digital surface model (DSM), which includes all surface features. After describing the method, we then use it to document changes in a slow-moving landslide near Taylor Mill in northern Kentucky over a 14 year period using a combination of (1) regional airborne lidar coverage acquired by government agencies before our study began and (2) site-specific structure-from-motion (SfM) and uncrewed aerial vehicle (UAV)-lidar surveys acquired as part of our work. To calculate thresholds of detectable change for each pair of DEMs, we used the statistics of noise for the difference map in areas outside the landslide, where no elevation change is reasonably inferred to have occurred. Elevation changes within the landslide are then measured relative to these no-change areas. As we explain, this is a pragmatic approach useful in many cases of practical interest for which control points are limited or non-existent. We found that combinations of surveys yield a threshold of detection of ±0.05 to ±0.20 m, depending on the combination, and that vertical changes in the landslide above the threshold of detection were found for every combination of surveys, even as short as SfM surveys performed 2 weeks apart. This approach is broadly applicable to other study areas with topographic and bathymetric data acquired using various technologies and of differing quality. Thus, this technique provides a method to utilize older survey datasets that may be regional scale and lower resolution in combination with site-specific high-resolution datasets to detect real elevation change over time.

We conclude that whereas the most significant movement of the Taylor Mill landslide, which formed on a slope altered by human activity, was triggered by record precipitation, the landslide has continued to move and expand laterally onto a previously unaltered slope. With climate change increasing precipitation regionally (EPA 2016), it is expected that landslide activity in the region will only increase, making the monitoring of slow-moving landslides a critical part of landslide mitigation.

The Cincinnati, Ohio and northern Kentucky region (Fig. 1) is plagued by slow-moving landslides that regularly threaten infrastructure (Mirus et al. 2017). For example, several 2019 landslides occurring along a major road cost over $17 USD million and took 2 years to remediate (City of Cincinnati Transportation and Engineering 2019; Knight 2021). Landslides typical of the region include slow-moving debris slides, generally <2 m thick, that form in the colluvium, unconsolidated material weathered from the underlying bedrock, which covers the steep slopes along the Ohio River and tributary valleys (Varnes 1978; Fleming and Johnson 1994); and deeper, slow-moving slumps that occur on flatter slopes in thick colluvium, glacial deposits or anthropogenic fill (Baum and Johnson 1996). The colluvium consists of weathered rock fragments ranging in size from granules to tabular limestone boulders in a clayey matrix. The colluvium most susceptible to landslides is derived from horizontally bedded shale and limestone of the Ordovician Kope Formation (Fig. 1), which is composed primarily of weak illitic shale that slakes easily when exposed to water (Koralegedara and Maynard 2017). The landslide chosen for this study, referred to as the Taylor Mill landslide (39.034234, −84.512587), is a translational debris slide that has been active since at least 2003 and has affected a slope and roadway leading to an apartment complex (Fig. 2). Attempts to mitigate the slide have included regrading the slope and rebuilding part of the roadway (Fig. 3). In 2012 the debris slide measured c. 45 m wide and 70 m long and by 2021 it had more than doubled its width. If we assume a constant thickness of 1.5 m, the volume of the landslide was initially c. 4500 m3 and has expanded to 11 500 m3.

Change detection

Aerial and satellite platforms have made it feasible to detect changes in slow-moving landslides with greater spatial and temporal resolution than with traditional methods such as field-based mapping or repeated surveys of benchmarks (Turner et al. 2015; Schulz et al. 2017; Okyay et al. 2019; Lacroix et al. 2020b). In particular, UAV-based lidar and SfM photogrammetry have made data acquisition less expensive and easily repeatable over short time periods (Jaboyedoff et al. 2012). Because slow-moving landslides may be active for many decades, older elevation datasets (e.g. topographic maps, Shuttle Radar Topography Mission (SRTM) or regional lidar) can also provide valuable past landslide information, albeit at lower spatial and temporal resolutions than site-specific SfM or lidar surveys.

A critical challenge in the study of elevation change detection, including slow-moving landslides, is addressing signal-to-noise ratios (Wheaton et al. 2009; Schaffrath et al. 2015; James et al. 2017). Frequently monitored landslides often register changes so small that the elevation changes may not be greater than the survey errors, referred to as noise. There are many potential sources of DEM noise, such as the incomplete density of observations, processing errors, measurement errors and errors introduced by the interpolation of point cloud data (Wechsler 2007), all of which will propagate through derivative maps and calculations (Holmes et al. 2000; Haneberg 2006, 2008). Survey noise can be quantified for SfM surveys (Clapuyt et al. 2016; Goetz et al. 2018) and elevation data acquired from various spaceborne and airborne platforms (Haneberg 2006, 2008; Gonga-Saholiariliva et al. 2011). However, the noise for older elevation datasets may be unknown. Thus, combining surveys acquired with different technologies to detect decadal change in a landslide can be a challenge. Few studies have directly used multiple technologies acquired at different times to detect change (e.g. Warrick et al. 2019), and quantifying the minimum detectable change possible using combinations of technologies has rarely been explored (e.g. Warrick et al. 2017).

Although change detection is optimally based on comparing results obtained from remote methods with surveyed monument locations on the ground, it is not always practicably possible. For example, the study area may be inaccessible to surveyors, monuments may be vandalized or accidentally destroyed between surveys, or a study may use historical data that predate the installation of survey monuments. In our case, the first two lidar datasets were acquired 12 and 7 years, respectively, before our work at the Taylor Mill site began and no historical survey monument data were available. To compensate for the lack of survey monuments, we expand upon an approach based upon the statistics of inferred areas of no change proposed by Haneberg (2017), whereby we quantify the noise in areas outside the landslide where no elevation change has occurred, to calculate the threshold of detectable change. Elevation changes in the landslide above this threshold are measured relative to the no-change areas outside the landslide. If the no-change area is indeed moving, then the calculated changes in the landslide are not absolute. However, they will still be very useful for understanding the movement of the landslide.

We measured changes in the landslide surface between 2007 and 2021 using differences between pairs in a series of (1) DEMs from county-wide airborne lidar surveys conducted in 2007 and 2012, (2) SfM DEMs derived from photographs acquired by a UAV in 2019 and 2020 and (3) a lidar DEM created from a UAV-lidar point cloud acquired in 2021. Each DEM was adjusted to the 2012 lidar DEM by removing coherent noise in the form of bias (an elevation difference of the same magnitude across the map area) and tilt (a systematic change in elevation difference across the map area) as described in detail below. The 2012 lidar dataset was chosen instead of the 2007 lidar dataset because it has a denser point cloud and fewer obvious artefacts compared with the 2007 data. Once bias and tilt between datasets were removed, we calculated the threshold of detectable elevation change for each pair of DEMs using the statistics of errors in areas outside the landslide where no elevation change is thought to have occurred using two methods. The methods are summarized in Figure 4 and described in detail in the following paragraphs.

Data acquisition

The county-wide lidar data and aerial photographs used in our analysis were obtained during leaf-off conditions in the winter of 2007 and 2012 (Table 1). The 2007 lidar data covering northern Kentucky were obtained by the Northern Kentucky Planning Commission, and have an average point spacing of 1.6 m. The 2012 lidar data and photographs were acquired as a part of the Kentucky Aerial Photography and Elevation Data Acquisition program (KYAPED). The 2012 lidar data were collected at an average of 0.68 m point spacing or better and a required vertical accuracy of ±0.15 m or better. The 2012 aerial photographs have 0.15 m pixels. The 2012 and 2007 lidar data were originally delivered in Kentucky state plane coordinates in feet; we projected them from state plane coordinates to UTM Zone 16N coordinates in metres and converted the vertical units from feet to metres. The Kentucky state plane and UTM coordinate systems both use the same horizontal datum (NAD 1983) and ellipsoid (GRS 1980), thus errors introduced by this projection are negligible.

We acquired digital aerial photographs for the SfM surveys in 2019 and 2020 using Mavic Phantom 2 and Mavic Pro UAVs. Flights were planned with DroneDeploy software (Hinge et al. 2019) to ensure that the images had 75% front overlap and 70% side overlap.

Our UAV-lidar data were acquired using a Matrice 600 Pro UAV and Yellowscan Surveyor lidar system in December 2020. Flights were planned using UgCS software to follow the topography at an elevation of 40 m. The lidar flight extended farther to the south than the 2019 and 2020 SfM surveys to include a heavily vegetated area where the SfM surveys did not produce usable results.

Data processing

The photographs for the SfM survey were processed to produce point clouds using Agisoft Metashape (James et al. 2017), and georeferenced using 19 easily identifiable ground control points (GCPs) in areas outside the landslide where no change was inferred to have occurred. GCPs were selected using the coordinates of reference points such as sewer grates and light poles visible in the KYAPED 2012 aerial photographs and the corresponding elevations from the 2012 lidar-derived DEM. After image alignment, the total residual error for the GCPs is <0.02 m. Point clouds were classified in Metashape using the Classify Ground Points tool, using a maximum angle of 15°, maximum distance of 0.5 m and cell size of 50 m. The SfM point clouds have an average point spacing of 0.040–0.046 m (Table 1).

The UAV lidar data were post-processed using local Continuously Operating Reference Station positioning data to improve their spatial accuracy (Olsen et al. 2013). Point elevations were further improved using strip adjustment and classified into ground and non-ground points using Yellowscan CloudStation software. The point cloud has an average point spacing of 0.049 m.

Point clouds from the county-wide lidar, SfM and UAV-lidar were processed using ArcGIS to produce 0.1 m DEMs using natural neighbour interpolation, a technique used to construct a surface from irregularly distributed points, to fill voids (Sibson 1981). Although interpolating the two county-wide lidar DEMs to 0.1 m does not add information or increase the resolution of features, it does facilitate comparison with the more detailed SfM and UAV-lidar DEMs. Although geomorphological change can be quantified from digital elevation data using grid-based (DEM) or point cloud comparisons (e.g. Qin et al. 2016; Okyay et al. 2019), we used gridded DEMs in this research because calculating DEM differences is easily performed using map algebra within GIS software and has a long history of successful application in geomorphological change detection studies (e.g. Okyay et al. 2019). An area of interest (AOI) was used to define the processing extent so that all DEM grids would be aligned. Without this step, the DEM grid would originate at the southeasternmost point of each point cloud, and thus the grid for each dataset would be slightly different. All datasets are georeferenced in UTM Zone 16N, EPSG 26916.

Noise maps

Coherent noise in the forms of bias and tilt between each DEM and the 2012 DEM was visually assessed using noise maps. The 2012 DEM was chosen as the surface that all other DEMs were corrected to because it was better quality than the 2007 DEM, and also had high-quality aerial photographs taken at the same time, and thus allowed ground control points to be selected for processing the SfM data. Noise maps were symbolized to show apparent elevation differences only in the range of ±0.20 m to visualize just the noise in areas where no real elevation change between datasets is expected (no-change areas). The range of ±0.20 m was chosen because it produces a continuous display of noise values across the map area. The noise was then quantified by sampling the distribution of noise in the road above and below the landslide.

Correction for bias and tilt

Each of the DEMs we used required some correction of either bias (between lidar DEMs) or tilt (between SfM DEMs). Ideally, the mean elevation change in a no-change area should be zero; thus, a non-zero mean indicates a bias between the datasets. This bias was removed by adding or subtracting the calculated mean value from the DEM being corrected. To adjust the slight tilt observed in the SfM DEMs, the apparent elevation difference for each of a series of 30 points along the road above and below the landslide was used to create a correction surface using the Topo-to-Raster tool in ArcGIS, which uses an iterative finite difference interpolation technique. This correction surface was then subtracted from the DEM being adjusted to remove the tilt. The corrections for bias and tilt remove the coherent noise and the remaining random noise is then used to calculate the threshold of detectable change between the two datasets.

Threshold of elevation change detection

When two elevation surveys are noisy and are then combined to calculate an elevation change, the errors are larger than the sum of the errors in each individual survey, which is often referred to as the propagation of errors (Birge 1939). Previous researchers have used probabilistic geomorphological change detection thresholds based upon the propagation of elevation errors from each of the two DEMs being compared (Brasington et al. 2000; Lane et al. 2003; Schaffrath et al. 2015). However, there may be cases of practical interest in which the individual DEM elevation errors are not available. For example, it is generally not feasible to directly determine DEM errors in deep-water seafloor change detection studies based upon repeat multibeam echosounder surveys (e.g. Haneberg 2018). Likewise, it is not unusual to encounter situations such as the one we describe in this paper, in which some of the DEMs were produced before our research began. Even if quality assurance statistics are available for previous topographic surveys, those data are typically collected in unobstructed, smooth and flat areas that are not representative of heavily vegetated, rough and steep landslide terrain. Such quality assurance error statistics can be misleading because DEM errors in areas relevant to landslide studies can easily be an order of magnitude larger than those collected for quality assurance purposes (e.g. Haneberg 2006, 2008). To account for these difficulties, we use two methods that combine error propagation theory with geomorphologically reasonable assumptions about areas in which no change is expected to have occurred. The elevation changes in the landslide that we report are relative to the no-change areas, which we assume have an elevation change of zero.

Both of our methods are based upon estimates of the propagated elevation difference error rather than the individual DEM elevation errors. This requires geomorphologically informed selection of one or more areas in which significant change can be reasonably inferred not to have occurred. We refer to these as no-change areas. Areas in the road above and below the landslide were used as no-change areas because there were no measurable elevation changes in these areas between 2007 and 2021 owing to either the nearby landslide or alterations to the road such as repaving. The propagated elevation difference error is estimated by calculating the mean and standard deviation of the DEM difference values within the no-change area(s) or, if the area is large, a subset of the values within the no-change area(s). In an ideal situation, the no-change data would be noise- or error-free; both the mean and standard deviation would be zero. Non-zero results for a no-change area are thus an estimate of the propagated elevation error. The mean, or bias, is removed to create a zero-mean no-change dataset. The remaining non-zero standard deviation, σΔz, is an estimate of the total propagated elevation error if the no-change inference is valid. The standard deviation of the propagated elevation errors is related to the individual DEM elevation errors by
σΔz2=σz,t12+σz,t22
(1)
where σz,t1 and σz,t2 are the standard deviations of the elevation errors in DEMs representing times t1 and t2, respectively (e.g. Hildebrand 1987; Brasington et al. 2000; Lane et al. 2003; Schaffrath et al. 2015). Next, we calculate the threshold value of detectable change using the propagated elevation difference error using two methods. Method I uses the statistics of errors from the difference map, whereas Method II uses back-calculated error estimates for each individual survey.

Method I: difference map errors

We assume that the DEM difference errors are normally distributed with a zero mean and heuristically adopt a threshold of ±2σΔz, so that the probability of a calculated elevation change being noise is <0.05 (Fig. 5a). This is essentially the same approach as taken by Brasington et al. (2000), Lane et al. (2003) and Schaffrath et al. (2015) except that equation (1) allows us to use the no-change area estimate of propagated error rather than the individual DEM errors. We also round the multiplier up from 1.96 to two for convenience.

Method II: individual DEM errors

We calculate a threshold value that must be exceeded for the overlap between the two DEM error distributions to remain below a specified significance level, α, as shown schematically in Figure 5b. If both component error distributions are known, α can be calculated numerically for any kind of distribution. For the work we describe in this paper, neither distribution is known. If it is reasonable to assume that both distributions are normal and have the same standard deviations, σz=σz,t1=σz,t2, either of the cumulative distribution functions can then evaluated at a value of one-half the threshold and the result doubled to obtain
Δzcrit=23/2σzerfc1(α)
(2)
where Δzcrit is the elevation difference threshold that must be exceeded to ensure the overlap between the two distributions is less than α, σz is the standard deviation of the individual DEM errors, and erfc−1 is the inverse complementary error function. For a normal distribution of errors with a mean of zero and a standard deviation of 1/2, the error function erf z=22π0zet2dt gives the probability that the error lies within ±x; the inverse complimentary function erfc−1 uses the probability to find ±x. Assuming that both DEM error standard deviations are equal, equation (1) can be rearranged to yield σz=σΔz/2 and equation (2) can be rewritten in terms of the difference map standard deviation rather than the individual DEM standard deviation:
Δzcrit=2σΔzerfc1(α).
(3)
For the significance level α = 0.05, or 95% confidence level, which is commonly used in many scientific studies, Δzcrit=2.77σΔz. Once the threshold values for each DEM combination are calculated, threshold maps are symbolized with a neutral colour representing values under that threshold.

Magnitude of vertical change

We use the magnitude of vertical change between the corrected DEMs integrated over a specified sample in the centre of the landslide area to quantify the amount of vertical surface deformation between each pair of DEMs. The magnitude of vertical change is the sum of the absolute value of the maximum change in each 0.1 m by 0.1 m cell in the sample area. We did not conduct a volume analysis on the landslide because unknown volumes of material have been added to and removed from the slide. An unknown volume of fill was placed on the top of the slope during road reconstruction between 2007 and 2012. As the landslide has repeatedly run out over the roadway, the toe has been excavated multiple times, removing an unknown volume of material.

DEM corrections

Each dataset required correction to vertically align it with the 2012 DEM (Fig. 6, Table 2). The 2007 county-wide lidar DEM and the 2021 UAV-lidar DEM were corrected for a bias of 0.057 and 0.058 m, respectively. All three SfM DEMs were corrected for tilt of up to ±0.21 m across the map area. Figure 6 shows examples of noise maps and noise distributions in no-change areas outside the landslide. Figure 6a and b shows an example of the noise map before and after correction for bias between the 2012 county-wide lidar DEM and the 2021 UAV-lidar DEM. Figure 6c and d shows an example of the noise map and noise distribution before and after correction for the tilt between the 2019 SfM DEM and the 2012 lidar DEM.

Threshold of detectable change

Table 2 shows thresholds of detectable change prior to and after corrections for bias and tilt for all DEM combinations. In 10 of the 15 combinations, the threshold of detectable change was reduced by up to 0.20 m after corrections; in two combinations, the threshold of detectable change increased by <0.01 m, and in three combinations, there was no change in the threshold value. After corrections, the threshold of detectable change ranges from ±0.05 to ±0.20 m, with the largest thresholds resulting from combinations that include the 2007 county-wide lidar DEM.

The threshold of detectable change calculated using Method I is always inherently smaller than that for Method II. For our datasets, the difference between the two methods ranges between 0.02 and 0.09 m, depending on the combination (Table 2). Threshold maps do not change underlying data, just the way it is symbolized, with a neutral colour representing elevation change under the threshold value. Figure 7 shows an example of the threshold maps resulting from Method I (Fig. 7a) and Method II (Fig. 7b) applied to the 2019 and 2020 SfM DEMs. Method I results in a calculated threshold of 0.05 m, and Method II results in a threshold of 0.07 m.

Vertical changes in the landslide surface

Real elevation change in the landslide is detected in every combination of surveys, including those flown just 2 weeks apart. Elevation changes in the landslide are measured relative to the no-change areas, which have an assumed elevation change of zero. Figure 8 shows a sequence of maps showing elevation changes in the landslide between 2007 and 2021 using a gradational scale. Elevation profiles that show the progressive elevation change between 2007 and 2021 in the landslide are shown in Figure 9. Observations of elevation change in the landslide are summarized as follows.

  • Between 2007 and 2012 (Figs 8a and 9a), a landslide c. 45 m wide and 70 m long developed in fill that had been placed on the slope. During this time, the head of the landslide dropped up to 2.2 m while the toe rose up to 1.3 m.

  • Between 2012 and 2019 (Figs 8b and 9b), most of the landslide continued to change in elevation. The top half of the slide shows some bands of elevation gain and loss that parallel the slope contours, whereas the toe mostly gains in elevation up to 1.2 m. In addition, a new area to the south lowered in elevation up to 1.2 m,although the lateral extent is obscured by vegetation and the topography could not be determined from the SfM data.

  • Between 2019 and 2020 (Figs 8c and 9c), the top half of the slide again shows parallel bands of small elevation gains and losses. Larger elevation changes are seen in the lower half of the slide, including the excavation of the toe near the road, which lowered the elevation there up to 1.8 m, whereas the northern part of the toe gained elevation of up to 0.9 m.

  • Between 2020 and 2021 (Figs 8d and 9d), elevation changes of up to 1.1 m are seen in the lower half of the landslide, whereas the upper half of the landslide appears to have little to no vertical change.

  • Difference maps showing the elevation change between lidar DEMs in 2012 and 2021 (Fig. 8e) allow a view of the ground beneath the vegetation and show that the landslide has extended to the south by at least 60 m along the lower portion of the slope. Gullies can be observed above, below and across the body of the landslide. Other areas of smaller elevation changes are seen to the north of the slide boundary, and in the mid-slope area above the southern part of the slide.

The sum of the positive elevation gains and the absolute value of negative elevation losses inside the landslide sample area for each combination of DEMs is shown in Figure 10, along with annual rainfall between 2005 and 2021. The rate of elevation change is greatest between 2007 and 2012, which coincides with the annual rainfall in 2011 of 1.86 m, the greatest annual rainfall since precipitation records have been collected starting in 1871. The rate of elevation change is also greater between 2019 and 2021. When positive and negative elevation changes are viewed individually (Fig. 10c and d), the rate of elevation gain is similar in all combinations, whereas the rate of elevation loss is greatest between 2007 and 2012, and between 2019 and 2021.

Threshold of detection

Establishing thresholds of detectable change is important for two reasons. First, it provides a way to utilize older or noisy low-resolution datasets in combination with more recent higher-resolution datasets to detect real change and offer information about the behaviour of slow-moving landslides over many years to decades. In this study, without the 2007 data, one might recognize a landslide from the topographic signature in the 2012 data, but we would not have information about the timing or magnitude of the change in the landslide. The methods of calculating the threshold presented here are broadly applicable to a wide range of topographic and bathymetric change detection problems using combinations of DEMs from different sources and for which limited control points may be available. Using the statistics of noise in no-change areas from the difference map, one can simply use Method I to calculate the threshold of detectable change. Method II, which always produces a higher threshold than Method I, can also be used; however, Method II uses the assumption that the standard deviation of noise is equal for both datasets. If the older data have a markedly different standard deviation of errors, equation (1) can be used to calculate the propagated error, which would require knowledge or estimation of the standard deviations, or at least the ratio of standard deviations for each dataset. This propagated error can then be used in equation (3) to calculate the threshold of detectable change.

Second, calculating a threshold of change can allow for more reliable and cost-effective monitoring of continuing changes in slow-moving landslides. In this study, annual changes to the landslide are easily detectable with either SfM or lidar, and the slight corrections we demonstrated here to reduce the threshold of detection may not even be necessary to accomplish this. However, in the case where critical infrastructure may be damaged by a slow-moving landslide, monitoring for small changes may be a crucial mitigation strategy; thus, a technique to quantify and minimize the threshold of detectable change, including making the small corrections for tilt and or bias, is important.

An example of a small change in the landslide that is detectable only after making corrections and using the resulting threshold is provided by the combination of surveys collected just 2 weeks apart on 11 March 2019 and 26 March 2019. Before corrections to the individual DEMs, the threshold of detectable change was 0.33 m (Method II), in which case there was no detectable change in the slide beyond the excavation of the toe (Fig. 11a). Once corrections were made to each DEM, the threshold of detectable change was lowered to 0.13 m (Method II), and small changes within the landslide become apparent (Fig. 11b). A comparison with a combination of DEMs after an additional 9 months had passed shows that these small changes correlate with the pattern of elevation change over that longer time period (Fig. 11d).

The bias between lidar DEMs or tilt between SfM DEMs was easily identified using the noise maps. The pattern of tilt, in particular, would have been easily obscured without the noise maps. In this case, samples of noise in the road would show inconsistent distributions of errors across the map area; or if a single noise measurement was made of the area outside the slide, errors with a large standard deviation, and thus a larger threshold of detection, would be calculated.

When vertical changes in the landslide are large compared with the noise, they are easy to distinguish even without reliance on a formally calculated threshold. When the changes are small, however, minimizing the threshold of detection is critical for distinguishing real change.

Landslide deformation

The series of difference maps between 2007 and 2021 show that the Taylor Mill debris slide developed on a slope that had been altered by road construction, regrading and the addition of fill. There is some evidence that a slow-moving landslide may have already existed at this location, but the most significant landslide movement occurred by 2012, probably owing to record precipitation in 2011. Since 2012 the lower half of the landslide has continued to slowly deform and move out into the roadway, where it is periodically excavated, whereas the scarp region displays only minor elevation changes. The most significant change since 2012 is the major lateral expansion of the landslide onto the natural slope, which had not been previously altered or shown any signs of movement. This demonstrates that a landslide that arguably was staged by human activities and triggered by record rainfall can lead to the destabilization of adjacent undisturbed slopes. This is important for the Cincinnati and northern Kentucky region because there are numerous other slow-moving landslides in the area that are similarly sensitive to changes in precipitation (Sparling 2019). Regional precipitation has increased 5–10% over the past 50 years owing to climate change (EPA 2016), and an increase in landslide activity is expected (Leshchinsky et al. 2017; Mirus et al. 2017; Coe et al. 2018). Furthermore, these landslides can potentially expand laterally onto previously undisturbed slopes, as seen at in this study, causing additional damage to property and infrastructure.

Imagery accessed via Google Earth prior to the first lidar survey of 2007 indicates that a landslide may have already existed on the slope, as imagery from 2003 and 2005 shows some bare patches in the otherwise vegetated slope, and the vegetation was cleared and slope regraded by June 2006 (Fig. 3). The county-wide lidar DEM and imagery from 2007 shows this regraded slope, and there is no evidence of the landslide at that time. Between May 2007 and October 2008 fill was added to the upper slope, and between October 2008 and July 2010 the upper curve of the driveway was reconstructed, and fill was again added to the upper slope (Fig. 3).

The major elevation changes within the landslide between 2007 and 2012 (Figs 8a and 9a) were probably triggered by the record rainfall of 2011 (Fig. 10b). The magnitude of elevation loss in the scarp area (up to 2.2 m) is much larger than the magnitude of positive elevation gain in the lower half of the slope (up to 1.3 m), as can be seen in the elevation profiles and the magnitude of vertical change (Figs 9a and 10). The topographic profiles, mapped elevation changes and magnitude of vertical change between 2007 and 2012 all indicate that more material was lost at the scarp than gained at the toe. Because the landslide toe extended into the road, this material was excavated and moved off-site.

Between 2012 and 2019 elevation changes continued to occur throughout the landslide, but the changes had a smaller magnitude (maximum of 1.2 m elevation change) than those that had occurred between 2007 and 2012 (Figs 8b and 9b). Parallel bands of elevation gain and loss that roughly parallel contours indicate either internal deformation of the slide material as it moves over an irregular slip surface or small internal slumps. Between 2019 and 2021 elevation changes have been greatest in the lower half of the slide, which has continued to advance into the driveway, and where debris has been periodically excavated (Figs 8c, d and 9c, d). Elevation changes in the upper half of the slide appear to be confined to small slides or erosion. Between 2020 and 2021, the upper half of the slide appears quiescent (Figs 8d and 9d).

By 2019 the landslide had also expanded laterally to the south. Although the full extent of this expansion was obscured by vegetation in the SfM surveys in 2019 and 2020, the lidar DEM in 2021 revealed that the landslide had laterally expanded by at least 60 m, thereby more than doubling its width (Fig. 8e). The landslide is possibly expanding to the north as well. Some of this expansion appears to involve the fill placed on the slope before the 2012 lidar survey but has extended significantly beyond that. The 2021 lidar DEM (Fig. 8e) also shows that several gullies had developed on the slope, which were not observed in the 2012 lidar, indicating the contribution of surface water to the landslide. Thus, the mobilization of the landslide by 2012 was followed by continued deformation of the toe, and a significant lateral expansion into a formerly undisturbed natural slope.

Landslide signature

This pattern of elevation change seen in the Taylor Mill landslide is characteristic of a debris slide, which ideally has an even loss of elevation in the scarp area, internal deformation in the body of the slide and elevation gain in the toe (Fig. 12). Interbedded shale and limestone produce an uneven slip surface, which helps generate internal deformation of the debris slide material (Fleming and Johnson 1994). A slump, in contrast, would have the greatest elevation drop at the scarp, which would progressively diminish to the axis of rotation of the slide where the change should be zero, and then a gradual increase in elevation towards the toe (Fig. 12). Thus, we observe that the mapped elevation change signature reveals the nature of the slide.

Because the landslide toe where it has run out over the road has been excavated periodically, we do not observe as much elevation gain in the toe as might be expected for a pristine debris slide. Typical debris-slide thicknesses are up to 2 m for landslides in the region (Fleming and Johnson 1994; Baum and Johnson 1996), and the slightly higher elevation loss values observed in the scarp may be due to the regrading of the slope in 2006, and the additional fill material that was added to this slope between 2007 and 2010.

Although differencing two DEMs produces a map of vertical change, this does not necessarily imply that the ground deformed only vertically. Elevation change in a landslide could be produced by vertical change but also by horizontal movement, or a combination of the two. In a translational debris slide, for example, there could be significant horizontal movement of the slide down the slope, whereas the vertical changes may not appear significant (Fig. 12). Continuing work on landslides in the Cincinnati and northern Kentucky region will address the horizontal component of movement, in addition to the vertical changes reported here.

Thermal expansion and contraction

The amount of expansion and contraction for diurnal or seasonal temperature changes will depend on the soil characteristics, moisture conditions and the temperature conditions. A recent study on a loess slope in China found that the soil expanded and contracted with an amplitude of about 1 mm over the course of a year owing to temperature variations (Lan et al. 2021). In a laboratory study of bentonite clay under thermal loading, the heating and cooling resulted in expansion and contraction with a volumetric strain of c. 1.0% (Tang et al. 2008). If we consider a hypothetical 2 m thick volume of bentonite that is allowed to expand only vertically as a result of thermal loading, an elevation change of c. 0.02 m would result. The clayey soils in the northern Kentucky area probably fall somewhere between these two examples, thus expansion and contraction values, although non-zero, would lie below the threshold values calculated for this study, the smallest of which is 0.05 m.

Influence of vegetation

The influence of vegetation on SfM DEMs is well documented (Cook 2017; Zekkos et al. 2018) and is significant when comparing elevation maps over longer time spans when changes to vegetation can be expected. However, we found that in SfM DEMs conducted close in time and prior to the growing season, very little change to vegetation has occurred. Therefore, any vertical change to the landslide will change the vegetation as well, and vertical changes in the vegetation reflect real vertical change in the landslide. The combination of DEMs over a 2 week period was able to detect small elevation changes in the body of the slide, and the locations and magnitudes of these changes were confirmed in the survey of the following year (Fig. 11).

It is sometimes assumed that vegetation is completely removed in bare earth DEMs derived from lidar, when in fact vegetation may still influence the DEM. For example, elevation changes between two lidar DEMs are seen in the mowed area across the road from the landslide, which reflects the differing grass heights during these two surveys (Fig. 6b).

Spatial variability of uncertainty

There is spatial variability in uncertainty related to terrain characteristics such as roughness (Podobnikar 2016), slope (Xiong et al. 2018), point density and vegetation (Clapuyt et al. 2016) and combinations of terrain characteristics (Carlisle 2005). In this study we used errors in the apparently stable and smooth roadways outside the landslide. If no road or other smooth surface was available, we expect that thresholds of detectable change would be higher. For example, for the combination of the 2012 county-wide lidar DEM and the SfM DEM of March 11, 2019, the road areas had a distribution of noise of 0.034 ± 0.033 m, and threshold values of 0.06 m (Method I) and 0.09 m (Method II). In contrast, the distribution of noise in the grassy area below the landslide was −0.074 ± 0.072 m, which produced threshold values of ±0.14 m (Method I) and ±0.20 m (Method II).

Slow-moving landslides are chronically destructive and can permanently damage property and infrastructure. The purpose of this study was to better understand decadal-scale spatial and temporal variations in a slow-moving debris slide in northern Kentucky over 14 years. To accomplish this, we needed to integrate existing county-wide lidar data acquired before our study began, along with site-specific SfM and UAV-lidar surveys. Because the technology and quality differ between surveys, the challenge was to devise a method to quantify survey noise so that a threshold of detectable change could be calculated. To reduce the threshold of detectable change, bias between lidar DEMs and tilt between SfM DEMs was first corrected to produce a zero mean elevation difference in areas outside the landslide where no change is inferred to have occurred. The threshold of detectable change was then calculated from the remaining random noise using two methods, each of which used the propagated elevation difference errors for the DEM combination. Method I used the errors from the difference map and Method II used back-calculated estimates of the individual DEM errors. The thresholds of detectable change range from ±0.05 to ±0.20 m, depending on the DEM combination and method used, with Method II producing a larger threshold value.

The series of difference maps between 2007 and 2021 show that the landslide developed on a slope that had been altered by road construction, regrading and the addition of fill. The greatest change in the landslide occurred between 2007 and 2012 and was probably triggered by the record rainfall of 2011. Since 2012 the lower half of the slide has continued to slowly deform whereas the upper half of the slide has been generally quiescent. The most significant change since 2012 is that the landslide has expanded laterally by at least 60 m into an unaltered slope that had previously shown no signs of movement, demonstrating that a landslide that was staged by human activity can lead to the destabilization of adjacent slopes.

The region has numerous other slow-moving landslides, and with precipitation increasing because of climate change, landslide activity is projected to increase. Thus, monitoring these slow-moving landslides is a critical part of landslide mitigation. We found that real change occurred in the landslide in all DEM combinations, including in SfM DEMs separated in time by just 2 weeks. Thus, in addition to providing a way to utilize older or noisy low-resolution datasets to document the behaviour of a landslide over many years to decades, the methods of calculating a threshold of detectable change presented here can also provide a reliable method of monitoring continuing changes in a slow-moving landslide. This technique for calculating a threshold of detectable change is widely applicable to other change detection studies where various survey technologies have been used to capture elevation data, including older regional or low-resolution surveys that might predate a particular investigation.

We thank J. Dortch, R. Thigpen, Y. Zhu, R. Olsen, and R. Noble-Varney for their contributions to the preparation of this paper.

SEJ: conceptualization (equal), data curation (lead), formal analysis (equal), investigation (lead), methodology (equal), project administration (lead), software (lead), visualization (lead), writing – original draft (lead), writing – review & editing (equal); WCH: conceptualization (equal), formal analysis (equal), methodology (equal), validation (lead), writing – original draft (supporting), writing – review & editing (equal); LSB: conceptualization (supporting), visualization (supporting), writing – review & editing (equal); MMC: investigation (equal), writing – review & editing (equal)

This research received no specific grant from any funding agency in the public, commercial or not-for-profit sectors.

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

The 2007 county lidar data can be found at https://apps.nationalmap.gov/downloader/. The 2012 county lidar can be found at https://kygeonet.ky.gov. Site-specific SfM and lidar data are available upon request.

This is an Open Access article distributed under the terms of the Creative Commons Attribution 4.0 License (http://creativecommons.org/licenses/by/4.0/)