Articles | Volume 19, issue 3
https://doi.org/10.5194/nhess-19-679-2019
https://doi.org/10.5194/nhess-19-679-2019
Research article
 | 
29 Mar 2019
Research article |  | 29 Mar 2019

Changes in ground deformation prior to and following a large urban landslide in La Paz, Bolivia, revealed by advanced InSAR

Nicholas J. Roberts, Bernhard T. Rabus, John J. Clague, Reginald L. Hermanns, Marco-Antonio Guzmán, and Estela Minaya
Abstract

We characterize and compare creep preceding and following the complex 2011 Pampahasi landslide (∼40 Mm3±50 %) in the city of La Paz, Bolivia, using spaceborne radar interferometry (InSAR) that combines displacement records from both distributed and point scatterers. The failure remobilized deposits of an ancient complex landslide in weakly cemented, predominantly fine-grained sediments and affected ∼1.5 km2 of suburban development. During the 30 months preceding failure, about half of the toe area was creeping at 3–8 cm a−1 and localized parts of the scarp area showed displacements of up to 14 cm a−1. Changes in deformation in the 10 months following the landslide demonstrate an increase in slope activity and indicate that stress redistribution resulting from the discrete failure decreased stability of parts of the slope. During that period, most of the landslide toe and areas near the head scarp accelerated, respectively, to 4–14 and 14 cm a−1. The extent of deformation increased to cover most, or probably all, of the 2011 landslide as well as adjacent parts of the slope and plateau above. The InSAR-measured displacement patterns, supplemented by field observations and optical satellite images, reveal complex slope activity; kinematically complex, steady-state creep along pre-existing sliding surfaces accelerated in response to heavy rainfall, after which slightly faster and expanded steady creeping was re-established. This case study demonstrates that high-quality ground-surface motion fields derived using spaceborne InSAR can help to characterize creep mechanisms, quantify spatial and temporal patterns of slope activity, and identify isolated small-scale instabilities; such details are especially useful where knowledge of landslide extent and activity is limited. Characterizing slope activity before, during, and after the 2011 Pampahasi landslide is particularly important for understanding landslide hazard in La Paz, half of which is underlain by similar large paleolandslides.

Dates
1 Introduction

Creep – steady, imperceptibly slow movement under the influence of gravity – precedes many large landslides in both bedrock (e.g. Chigira and Kiho, 1994; Kilburn and Petley, 2003) and soil (e.g. Kalaugher et al., 2000; Petley et al., 2002). Its temporal patterns can improve understanding of slope-deformation processes (Petley et al., 2002) and the possible timing of impending catastrophic failure (Voight, 1989; Crosta and Agliardi, 2003; Federico et al., 2015). Its spatial patterns indicate the style (Gischig et al., 2011; Schlögel et al., 2016) and approximate magnitude of instability. Consequently, characterization of temporal and spatial patterns of slow deformation is a focus of many landslide hazard assessments, monitoring systems, and warning programs (Eberhardt, 2012; Baron and Supper, 2013; Hermanns et al., 2013a). In many cases, however, pre-failure creep is not quantified, because either it is not recognized prior to catastrophic failure or failure carries a relatively low risk to people and property.

Creep is rarely monitored following a landslide, except at sites where renewed catastrophic failure is expected (e.g. Gischig et al., 2011). Limited post-event monitoring sometimes stems from an assumption that stress release during catastrophic failure enhances slope stability. In addition, instrumentation that may have been present on or within the slope is typically destroyed during failure (e.g. Crosta and Agliardi, 2003). Consequently, activity soon after failure is generally unclear and can rarely be compared with pre-failure deformation. Such records are, however, of utmost importance in light of observed spatiotemporal clustering of failure events (e.g. Hermanns et al., 2006, and references therein; Crosta et al., 2017; Hilger et al., 2018).

Spaceborne radar interferometry (InSAR) has long been recognized for its promise in characterizing the spatiotemporal evolution of slopes (e.g. Petley et al., 2002; Colesanti and Wasowski, 2006; Wasowski and Bovenga, 2014). Continued technological advancements and data availability offer the potential for further improvements to the characterization of slope behaviour (Wasowski and Bovenga, 2014), including prediction of future failure (Moretto et al., 2017). Improvement of radar systems is increasing ground resolution and phase quality, and the establishment of satellite constellations is improving temporal resolution and interferometric coherence. Due to such improvements, pre-failure accelerated creep has recently been detected and measured, for both natural (Intrieri et al., 2018; Handwerger et al., 2019) and cut (Carlà et al., 2018) slopes. Such advancements similarly improve opportunities to analyze post-failure creep.

Our quantification of creep during the 30 months preceding and 10 months following the largest modern landslide in the city of La Paz, Bolivia (2011 Pampahasi landslide, ∼40 Mm3), using InSAR demonstrates complex slope behaviour. Contrary to the sometimes-invoked view that a major failure reduces large-scale instability, the affected slope shows enhanced post-failure activity, highlighting the complexity of slope stability in La Paz. The Pampahasi case study illustrates the necessity for land-use planning views that better align with the complexity and commonly recurrent nature of large landslides. The InSAR technique we employ is optimized to maximize displacement record density and distinguish abrupt changes in motion, which together improve the characterization of spatially complex deformation patterns typical of landslides. The resulting displacement fields provide nearly spatially continuous quantification of motion over large portions of the landslide and surrounding terrain that enable discrimination of localized motion and aid interpretation of failure mechanisms.

2 Setting and material properties

2.1 La Paz region

A tributary of the Amazon River system penetrates the Cordillera Real ∼60 km southeast of La Paz and drains the eastern margin of the otherwise internally drained Altiplano plateau (∼4000 m a.s.l.) (Fig. 1a). The city of La Paz (population ∼0.8 million) is located within the upper part of this watershed where the valley system has incised into up to 800 m of weakly lithified Neogene and Pleistocene sediments shed from the Cordillera Real (Dobrovolny, 1962) (Fig. 1b). South of the city, incision has exposed folded and faulted Paleozoic to Mesozoic metasedimentary basement rocks forming the root of the range (Fig. 1b). The overlying Cenozoic fill sequence comprises over 600 m of predominantly fine-grained lacustrine and fluvial sediments of the weakly lithified Miocene–Pliocene La Paz Formation. It coarsens upward and toward the cordillera into more competent late Pliocene and Early Pleistocene glacial and proglacial sediments (Roberts et al., 2018) that are up to 400 m thick. The fill has been offset by numerous faults trending northwest and west (Dobrovlny, 1962; Lavenu, 1978; Roberts et al., 2018) (Fig. 1b).

https://www.nat-hazards-earth-syst-sci.net/19/679/2019/nhess-19-679-2019-f01

Figure 1Setting of La Paz, Bolivia. (a) Location and physiographic context. (b) Generalized geology after Dobrovolny (1962) and Anzoleaga et al. (1977), with faults from Dobrovolny (1962) and Lavenu (1978). Terrain is from the ASTER GDEM v2 produced by METI and NASA.

Download

The slopes and floors of the La Paz valley system are mantled by colluvium, including deposits of large (up to 50 km2; Dobrovolny, 1968) paleolandslides derived predominantly from the La Paz Formation (Dobrovolny, 1962; Anzoleaga et al., 1977) (Fig. 1b). Most of the paleolandslides date to the Late Pleistocene and Holocene (Dobrovolny, 1962; Hermanns et al., 2012). Historic instability is concentrated within paleolandslide deposits (e.g. Fig. S1 in the Supplement), particularly at their margins, likely reflecting low residual strengths and possibly differential stresses related to suspected ongoing creep. Failures < 1 Mm3 in size occur yearly and happen mainly during the rainy season (December–March) (O'Hare and Rivas, 2005; Roberts, 2016). In contrast, of the seven historic landslides larger than 1 Mm3, four happened during the dry season (April–November; Table S1 in the Supplement). Although general geomorphic and geologic characterization has been undertaken for some of the larger landslides (Dobrovolny, 1962; Anzoleaga et al., 1977; Quenta et al., 2007, 2008; Hermanns et al., 2012), detailed site investigations are lacking. Land-use decisions in La Paz have historically overlooked recurrent slope failures, with large landslide complexes being either unrecognized or repeatedly resettled after reactivations.

https://www.nat-hazards-earth-syst-sci.net/19/679/2019/nhess-19-679-2019-f02

Figure 2Ground surface conditions in the Pampahasi area. (a) Surficial geology prior to the start of the twenty-first century; adapted from Anzoleaga et al. (1977). (b) Slope morphology and ground disturbance showing the extent and morphologic features of landslides that occurred during the period of radar acquisition (September 2008 to December 2011). Horizontal displacements during the 2011 Pampahasi landslide are based on comparison of pre-failure and post-failure WorldView-2 scenes. Contours are from NASA's Shuttle Radar Topographic Mission v3.0.

Download

2.2 Pampahasi area

Pampahasi Plateau is 2.2 km long, 0.2–0.7 km wide, and forms part of the interfluve between ríos Orkojahuira and Irpavi (Fig. 1b). Its margins are generally scalloped (Fig. 2) and in several places terminate in nearly vertical cliffs dropping 20–80 m to irregular, locally steep valley slopes below. The plateau surface dips slightly southward and lies at an elevation several hundred metres below the local southwest-dipping Altiplano and its erosional remnants (Fig. 3a).

The Río Irpavi valley is cut into well-stratified, fine-grained facies of the La Paz Formation, which comprises horizontal to sub-horizontal silt and sand beds interlensing with mainly granitic pebble–cobble gravel (Fig. 3b). The variable texture and limited (< 1 km) lateral continuity of these beds are typical of the unit, particularly its middle part (see Ahlfeld, 1945a; Dobrovolny, 1962; Bles et al., 1977). The age of the La Paz Formation beneath Pampahasi Plateau (Fig. 3b) is constrained by the 5.4 Ma Cota Cota Tuff (Lavenu et al., 1989; Servant et al., 1989), which crops out at ∼3520 m a.s.l. in the southern part of the city (Lavenu et al., 1989), and the 2.74 Ma Chijini Tuff (Roberts et al., 2017), which crops out along the west margin of the plateau (Dobrovolny, 1962) between 3800 and 3850 m a.s.l. The local Altiplano surface formed on top of the thick sequence of Pliocene–Pleistocene glacial sediments by 1.8 Ma or possibly as recently as 1.0 Ma (Roberts et al., 2018), attaining an elevation of approximately 4200 m a.s.l. in the vicinity of Pampahasi (Fig. 3a).

https://www.nat-hazards-earth-syst-sci.net/19/679/2019/nhess-19-679-2019-f03

Figure 3Stratigraphic sequence. (a) Cross-section oblique to the trend of the Cordillera Real through the fill underlying the Altiplano (location in Fig. 1b). (b) Lithostratigraphy beneath Pampahasi Plateau, including principle units involved in the Pampahasi paleolandslide and 2011 Pampahasi landslide. Modified from Anzoleaga et al. (1977) based on field mapping throughout the La Paz area.

Download

About 400 m of incision of the Altiplano surface by ancestral Río Irpavi after 1.8–1.0 Ma formed a south-sloping (Dobrovolny, 1962) valley-bottom unconformity onto which the Middle Pleistocene Pampahasi gravel (Fig. 3b) was deposited, as either fluvial channel deposits (Ahlfeld, 1945b) or an ancient landslide (Bles et al., 1977). This multilithic silty sandy pebble–cobble gravel is commonly more than 20 m thick (Dobrovolny, 1962) and may be up to 50 m thick along the plateau's southern edge where the unit overlies a locally steep erosional contact (Dobrovolny, 1955).

Many slopes descending from Pampahasi Plateau are mantled by ancient landslide deposits derived from the La Paz Formation (Dobrovolny, 1962; Anzoleaga et al., 1977). The deposits consist of largely intact brittle blocks of the La Paz Formation within a matrix of softer, ductile, and completely disaggregated La Paz Formation. Small (< 1 to 20 ha) failures evident in aerial photographs have occurred since at least the early twentieth century in many places within large (∼20–200 ha) prehistoric failures. The slope descending ∼380 m to Río Irpavi directly east of the southern half of Pampahasi Plateau comprises the undated Pampahasi paleolandslide (Fig. 2a).

Details of the prehistoric Pampahasi paleolandslide, which substantially predates historic records, must be inferred from its surface expression because no subsurface investigations have been conducted on this slope. Its morphology indicates a complex failure, but whether it represents one or multiple events is uncertain. Rotational displacement of the upper quarter of the slope transitions downward to eastward sliding and localized flow. The landslide displaced Río Ipravi eastward from its alignment farther upstream and downstream by at least 300 m. The substantial height of the head scarp (up to 80 m) suggests that the failure zone is many tens of metres below surface in the upper part of the landslide. The failure depth lower on the landslide is likely much shallower.

Surface incision of the Pampahasi paleolandslide by several gullies suggests that it is of substantial antiquity. The adjacent and similarly incised Villa Salomé paleolandslide is likely a separate landslide complex because its failure surface is higher and is separated from the Pampahasi paleolandslide along Río Jankopampa by a > 120 m thick sequence of intact La Paz Formation (Fig. 2a). In the area of the Pampahasi paleolandslide, the Río Irpavi channel is modified with check dams and concrete armouring. Before the 2011 failure, western tributaries of Río Irpavi draining the Pampahasi and Villa Salomé paleolandslides were little modified and had incised into weak loose material. The main channels draining the area (ríos Chujilluncani and Jankopampa, Fig. 2) have since been entombed in concrete culverts. In light of the geomorphic evidence of recurrent instability at the sites, Scanvic and Girault (1989) recommended that this area not be developed. However, the initially sparse development greatly expanded during the last decade of the twentieth century and first decade of the twenty-first century, resulting in the establishment of several large neighbourhoods.

2.3 Geomechanical properties

The different characteristics of the Pampahasi gravel and La Paz Formation give rise to contrasting material properties, which have been documented by Anzoleaga et al. (1977; summarized in Table S2). The Pampahasi gravel (50 %–75 % clasts, 19 %–50 % sand and silt, < 6 % clay) is permeable (k1×10-3 m s−1), has low to medium plasticity (plasticity index ≤10), and exhibits high internal friction (> 30), particularly when dry. Fine-grained zones are weaker but uncommon.

The properties of the La Paz Formation are variable due to the unit's heterogeneity. Coarse sand and gravel zones (35 %–75 % clasts, 13 %–65 % sand and silt, < 12 % clay) have permeability similar to that of the Pampahasi gravel but are commonly cemented by carbonate (≤15 %) or iron oxide. Fine-grained zones (≤10 % clasts, 40 %–90 % sand and silt, 10 %–40 % clay) are much weaker and, as the main facies, dominate the unit's behaviour. They are weakly cemented by carbonate (< 5 %) and locally abundant clay. The clay is largely non-expansive, but montmorillonite (smectite) is common in some beds and lenses. Due to low permeability (k<1×10-5 m s−1), water circulates mainly through fractures in interstratified gravel lenses or along joints. Plasticity is medium to high (plasticity index of 7–32), and shear strength ranges from 18 to 37. The unit is generally strong when dry and undisturbed (cohesion 0.2–0.7 MPa) and stands in near-vertical slopes. However, wet silty sediments have much lower cohesion (0.01–0.1 MPa) and shear strength. Smectite-rich clay lenses have greater cohesion (≥0.1 MPa) but high plasticity (plasticity index of 60–80) and low shear strength (≤14). Parts of the La Paz Formation that have failed in the past are even weaker due to loss of compaction and localized reduction of shear strength to residual levels (as low as 13 in silty and sandy zones and less than 8 in clay zones).

3 Methods

3.1 Failure event characterization

Our characterization of the 2011 landslide is based on field visits in the years before and after the event and on comparison of high-resolution optical satellite images acquired shortly before and after the event. Field observations by municipal staff and eyewitness accounts of residents provide details on the location, type, and magnitude of damage shortly before, during, and after the 2011 failure. We mapped features of the 2011 (26 February to 1 March) landslide from the first cloud-free, post-failure imagery (WorldView-2 acquired on 23 March 2011), and we quantified horizontal displacement vectors during the failure by comparing it with a pre-failure (5 January 2011) image from the same sensor, which we coregistered using points outside the landslide-affected area. Due to errors resulting from the lack of orthorectification and from elevation differences cause by the landslide, the estimated horizontal vectors are approximate (rounded to the nearest 5 m increment).

3.2 HDS-InSAR

We applied Homogeneous Distributed Scatterer InSAR (HDS-InSAR; Eppler and Rabus, 2011; Rabus et al., 2012) to a stack of 44 fine-beam-mode (5.1 m azimuth × 8.0 m ground-range pixel spacing) RADARSAT-2 ascending scenes (look direction: 36.3 incidence angle from vertical at scene centre, 76.0 aspect angle counter clockwise from north) acquired between September 2008 and December 2011 (Tables S5 and S6). We processed the scenes separately as a pre-failure stack (32 scenes) and a post-failure stack (12 scenes) to compare deformation patterns before and after the landslide and to reduce decorrelation caused by landslide-induced terrain changes. Although a range of displacement directions is possible, ground-surface displacement vectors due to mass movements are most probable in the downslope direction; we thus estimated true displacement rates from the angular deviation between downslope vectors and the satellite's line of sight (LOS).

The HDS approach combines the strengths of Persistent Scatterer InSAR (PS-InSAR; Ferretti et al., 2001) and Small Baseline Subset InSAR (SBAS-InSAR; Berardino et al., 2002) for considering, respectively, point and distributed targets. A continuously weighted, spatially adaptive filter defines spectrally similar pixel clusters (HDS neighbourhoods) within a rectangular search area based on similarity of their backscatter amplitude time series. The filter is applied to interferograms to generate differential phase and coherence for each HDS neighbourhood (Eppler and Rabus, 2011), which comprises either a single-pixel cluster (a persistent point target) or the weighted mean of a multi-pixel cluster (interpreted as the area covered by a homogenous distributed target). To more efficiently represent the high density of displacement records, we interpolated linear deformation maps from HDS to provide base images of long-term average LOS deformation over the entire area.

Like PS-InSAR, HDS-InSAR generates high-quality displacement time histories attributable to single pixels, each containing a dominating point target with near-linear displacement characteristics (Colesanti and Wasowski, 2006). Additionally, it provides time histories for multi-pixel clusters representing homogeneous distributed scatterers. Like SBAS-InSAR, HDS-InSAR characterizes strongly non-linear motion (Necsoiu et al., 2014) using the higher continuity coverage (Lauknes et al., 2010) of distributed targets. Unlike SBAS-InSAR, however, it does not significantly diminish phase quality or reduce spatial resolution. Two key differences set HDS-InSAR and some other advanced InSAR algorithms (e.g. Ferretti et al., 2011) apart from the PS-InSAR and SBAS-InSAR techniques: the use of adaptive filtering to preserve, as much as possible, spatial resolution while suppressing the noise from surface decorrelation of nearby incoherent pixels, and optimization to characterize spatially uncorrelated ground motion. These features make HDS-InSAR particularly well suited for characterizing landslides, which commonly display ground motion variability over short distances or across abrupt transitions. Although HDS-InSAR and SqueeSAR (Ferretti et al., 2011) share similar adaptive filtering methods, the latter uses wrapped phase triangulation to invert the interferometric network. In contrast, HDS-InSAR relies on prior unwrapping of the network interferograms, with potential unwrapping errors being corrected iteratively using the network redundancy a posteriori. Despite the differences in their algorithms, the final accuracy and spatial detail of both advanced InSAR methods are similar.

4 2011 Pampahasi landslide

The Pampahasi landslide1 occurred at the end of February 2011 following a week of variable precipitation (0–39.2 mm daily, Fig. 4), which included one of the wettest days on record (25 February, Table S3). The ∼1.5 km2 failure remobilized a large part of the Pampahasi paleolandslide deposit between the east side of Pampahasi Plateau and Río Irpavi (Fig. 2). In addition to previously failed material, it included small volumes of in situ weakly cemented La Paz Formation and, in the head-scarp area, the overlying uncemented Pampahasi gravel (Fig. 2a). None of the sediments involved were substantially weathered.

https://www.nat-hazards-earth-syst-sci.net/19/679/2019/nhess-19-679-2019-f04

Figure 4The precipitation record for the Laykacota meteorological station (see Fig. S1 for station location). (a) Comparison of monthly average precipitation (1945–2016) with monthly precipitation for the second half of 2010 and first half of 2011. (b) Daily precipitation and 7-day cumulative precipitation during February and early March 2011.

Download

Several locations in the upper part of the slope showed field evidence of slow ground deformation, largely visible as offset walls and road surfaces (Fig. 5, locations in Fig. 2a) in the years prior to catastrophic failure (Quenta and Calle, 2005; Hermanns et al., 2012). These features were recognized as evidence of an impending landslide (Quenta and Calle, 2005), although not of the magnitude of the 2011 failure. In response, the municipality of La Paz increased risk communication and installed concrete pillars to remediate what became the 2011 landslide head scarp (Hermanns et al., 2012). However, due to their localized and shallow nature, the stabilization efforts probably had little if any influence on the stability of the slope or behaviour of the subsequent failure.

https://www.nat-hazards-earth-syst-sci.net/19/679/2019/nhess-19-679-2019-f05

Figure 5Field evidence of pre-failure creep in what would become the head region of the 2011 Pampahasi landslide. (a) Small recent slumps (arrows) along Río Chujilluncani, ∼400 m downslope from the east margin of Pampahasi Plateau. (b) Minor cracking and centimetre-scale vertical offset of a paved road about 200 m east of, and 75 m below, the east margin of the plateau. (c) Unpaved road with ∼1 m displacement near the 2011 head scarp. (d) Vertical offset along a gully ∼75 m from the future south lateral margin of the 2011 landslide. The side of the gully nearer the middle of the 2011 landslide (right of the photo) has dropped ∼0.5 m relative to the adjacent slope, causing deformation of the gabions and curb in the foreground. (e) Cracks in an adobe wall (right foreground) and brick house (left middle ground) ∼50 m downslope of the east margin of the plateau, corresponding to the middle of the 2011 landslide head scarp. See Figs. 2a and 8a, b for locations of photographs. Photos taken by Reginald L. Hermanns in October 2005.

Download

The 2011 Pampahasi landslide occurred in several phases from the evening of 26 February to 1 March (Fig. 6, locations in Fig. 2b). Eyewitness accounts suggest that failure initiated in an area of earlier small slumps along Río Chujilluncani2 (Hermanns et al., 2012; Fig. 5a) ∼400 m downslope of the east margin of Pampahasi Plateau. Motion was fastest (up to several metres per second) and largely vertical in the head region, which failed shortly after, lasting 1 to 2 h and forming a ∼60 and 80 m high head scarp (Fig. 6a). Movement in the central zone and at the toe, which began soon after the initial formation of the head scarp, was slower (metres per minute to metres per day) and largely horizontal. The degree of ground disturbance and the magnitude of displacement generally decreased to the east and north (Fig. 2b). The northern part of the landslide's toe experienced relatively limited, but economically costly, ground disturbance (e.g. Fig. 6b–d) and moved as a more or less coherent mass ∼15 m southeastward (Fig. 2b). Deformation of the lower part of the landslide continued for days. A bridge across Río Irpavi (Fig. 6e) was shifted 10 m and uplifted on the morning of 27 February; uplift of several metres at the east end of the bridge suggests a rotational failure zone passing, at least locally, under Río Irpavi. Buildings near the south toe of the landslide (Fig. 6f) collapsed on 28 February and 1 March due to delayed relatively minor translation.

https://www.nat-hazards-earth-syst-sci.net/19/679/2019/nhess-19-679-2019-f06

Figure 6Impacts of the 2011 Pampahasi landslide. (a) Head scarp (80 m height) along the east margin of Pampahasi Plateau (note excavator and trucks in lower middle ground for scale). (b) Ground fissures in the northern part of the lower landslide mass. The wall in the background is ∼2 m high. (c) Back-tilted houses and apartment buildings on the northern part of the lower landslide mass (movement was from left to right). (d) Community building with collapsed wall on the northern part of the lower landslide mass. (e) Bridge across Río Irpavi that was deformed and shifted from its footings on 27 February near the centre of the landslide toe. (f) Buildings on the southern part of the toe of the landslide that collapsed on 28 February and 1 March (dates of collapse shown; photo by Marco-Antonio Guzmán, 2 March 2011). (g) Water service line ruptured by the landslide (upper pipe) and replacement water service line (lower pipe). See Figs. 2b and 8c, d for locations of photographs. Photos taken by Nicholas J. Roberts in May 2011 unless otherwise indicated.

Download

The surface morphology of the landslide likely reflects a complex failure zone with both rotational and translational components. Assuming a minimum average thickness of the landslide deposit of 15 m, its volume is at least 20 Mm3. However, the volume is more likely about 40 Mm3, given deeper failure in both the upper rotational part (40 m average depth over 32 ha) and the main predominantly translational body (25 m average depth over 110 ha).

Due to prompt evacuation of the head-scarp area and the dominantly moderate to low movement velocities farther downslope, no lives were lost during the landslide. However, about 1000 homes were destroyed and 6000 people were displaced. The rupture of a water line crossing the upper part of the landslide (Fig. 6g) left between 200 000 and 300 000 people in southern La Paz without potable water for several months (Hermanns et al., 2012; Aguilar, 2013). Continued episodic collapse of the steep head scarp in the years since the landslide has resulted in additional loss of homes and repeated expansion of an evacuation zone on Pampahasi Plateau. Large portions of the landslide complex have been recently resettled or are being prepared for reoccupation, with limited control of slope infiltration and runoff.

5 Deformation before and after failure

Due to the high spatial density of displacement records, particularly in high-coherence areas (Fig. 7b, e), interpolated linear deformation maps are preferable to HDS point data for interpreting the spatial variability of slope creep (see Fig. 7c, f). Furthermore, the interpolation weighting suppresses very localized HDS clusters that differ greatly from the average (Supplement, Sect. S4.7). Consequently, these small-scale variations, whether representing surficial movement or noise from uncorrected phase unwrapping errors, are insignificant to large-scale patterns described below. The displacement maps record deep, spatially regular slope movements as well as shallower, more variable movements, but their differentiation requires consideration of displacement patterns and may not always be clear.

https://www.nat-hazards-earth-syst-sci.net/19/679/2019/nhess-19-679-2019-f07

Figure 7Typical density and spacing of HDS-InSAR results for an area of mixed suburban, agricultural, and commercial development. Both the pre-failure stack (a, b, c) and post-failure stack (d, e, f) are shown. Panels for each stack depict satellite images showing land cover without InSAR data (a, d); point data representing HDS neighbourhoods (pixel neighbourhoods) (b, e), each recording a displacement time history; and linear displacement map interpolated from HDS points (c, f). Localized variability in HDS results is suppressed by linear deformation interpolation. Displacement rates are in the satellite line of sight. The deforming area is the southern boundary of the active part of the landslide toe. Base image is a 4 January 2011 QuickBird image viewed in GoogleEarth. Base images are from GoogleEarth (a, b, c 4 January 2011 and d, e, f 25 July 2011). See Fig. 8a, c for location.

Download

https://www.nat-hazards-earth-syst-sci.net/19/679/2019/nhess-19-679-2019-f08

Figure 8InSAR-measured line-of-sight ground displacements in the Pampahasi area. Average linear deformation rates during (a, b) the 30-month period before and (c, d) the 10-month period after the Pampahasi landslide. Compressed (−1 to 1 cm a−1; a and c) and expanded (−5 to 5 cm a−1; b and d) displacement scales emphasize, respectively, the spatial limits and spatial variability of slope activity. Dashed lines are the limits of landslides that occurred in 2009 and 2011 during the period of radar scene acquisition. Roman numerals show locations of displacement histories in Fig. 9. Letters show locations of photos taken, respectively, before (Fig. 5) and after (Fig. 6) the 2011 landslide. Base images are from GoogleEarth (a, b 4 January 2011 and c, d 25 July 2011).

Download

InSAR-measured ground deformation in the Pampahasi area is almost entirely restricted to mapped prehistoric landslide deposits, namely the Pampahasi and Villa Salomé paleolandslides (Fig. 2a), and the terrain directly behind their head scarps (Fig. 8). Their characteristic displacement histories are shown in Fig. 9 (locations i–x). Data gaps in each stack are locations of low coherence that mainly result from decorrelation related to remedial earthworks performed after a 10 ha landslide in 2009 (Fig. 8a, b) and the 2011 landslide (Fig. 8c, d) or from aliasing where displacement rates are greater than the detection threshold of RADARSAT-2 (∼2.8 cm LOS – equivalent to a two-way travel distance comparable to RADARSAT-2's 5.6 cm wavelength – over 24 days). Motion toward the satellite suggested by the data in the east-facing slopes (Fig. 8) is improbable, given that it implies uplift or motion upslope. Such issues are almost entirely restricted to the post-failure records (Fig. 8c), suggesting that the small size and limited redundancy of the interferometric network in the thinner stack has resulted in incomplete network correction of either phase unwrapping errors or atmospheric effects.

5.1 Pre-failure creep

Prior to the end of February 2011, about 180 ha of the lower part of the slope within the limits of the 2011 landslide was moving at rates as high as 14 cm a−1 (9 cm a−1 LOS; Fig. 8a, b). Motion was most widespread over the northern part of the toe of the Pampahasi paleolandslide (2–5 cm a−1 LOS, i in Fig. 8), where it extended into Río Jankopampa alluvium mantling the paleolandslide deposit. There it terminated abruptly ∼50 m north of the river (Fig. 8a) at a locally stable slope comprising undisturbed La Paz Formation (Fig. 2a). A 1 ha zone directly west of Río Irpavi and halfway along the toe of the paleolandslide (ii in Fig. 8) crept even faster (5–8 cm a−1 LOS, Fig. 8b). In contrast, activity on the southern part of the paleolandslide toe (iii in Fig. 8) was limited to a few small areas creeping at less than 0.5 cm a−1 LOS.

Pre-failure movement in the upper half of the Pampahasi paleolandslide was restricted to an area of no more than ∼250 m along the slope by ∼300 m downslope. Given the presence of stationary ground in some parts of the upper half of the landslide (north and south of B in Fig. 8a, b), creep may have been localized to a few smaller areas there. Minor isolated movement (< 1 cm a−1 LOS) is apparent along what became the 2011 landslide's upper lateral margins (e.g. iv in Fig. 8). Creep within the Pampahasi paleolandslide deposit beyond the limits of the 2011 failure was also localized and generally slower than 1 cm a−1 LOS (Fig. 8a). The largest of these areas involved 12 ha of paleolandslide material between the 2009 and 2011 landslides (v in Fig. 8) directly upslope (northwest) of the most active part of the slope toe and crept at 0.5 cm a−1 LOS. The most rapid pre-failure motion (9 cm a−1 LOS) was restricted to the upper portion of the future landslide where infrastructure damage had been documented in prior years (Fig. 5b, c). Aside from a small, slow-moving (0.5 cm a−1 LOS; vi in Fig. 8) area, no motion is evident along the future landslide head scarp, including sites of localized deformation observed in 2005 (Fig. 5d, e), or on the plateau behind it (vii–x in Fig. 8).

5.2 Post-failure creep

During the 10 months following the landslide, the area of movement on Pampahasi Plateau and the valley slope between it and Río Irpavi increased by about two-thirds to nearly 300 ha and included the entire area of the 2011 landslide where coherence was maintained (Fig. 8c). The maximum inferred downslope displacement rates following the landslide (14, 9 cm a−1 LOS; Fig. 8d) were similar to those before it but occurred in a region of the toe that was previously stable (iii in Fig. 8) and along the new head scarp (vi in Fig. 8) where creep had previously been slow (0.5 cm a−1 LOS). Displacement rates also increased in the northern part of the toe (ii in Fig. 8, from 2–5 to 3–6 cm a−1 LOS) and were more spatially variable than in the pre-failure period (Fig. 8).

Detection of movement is not possible over much of the middle and upper parts of the 2011 landslide due to decorrelation resulting from earthworks that continued for years after the event. However, because this zone is bordered on all sides by moving terrain (> 1 cm a−1 LOS, Fig. 8c), much of it also was likely creeping throughout the period of the post-failure stack. The rate and extent of creep increased in a zone measuring nearly 500 m by 500 m (including v) in the paleolandslide deposit between the 2009 and 2011 landslides. Given the limited transport, and thus bulking, across most of the area of the 2011 landslide (Fig. 2b) as well as the occurrence of abundant post-failure surface displacement beyond its limits (Fig. 8c), ground motions following the event should largely represent mass movements as opposed to soil settlement or compaction.

The most striking change in ground motion following the 2011 failure is the development of a new area of creep several hundred metres wide fringing the head region of the landslide (Fig. 8c). This zone extends up to nearly 400 m beyond the head scarp on Pampahasi Plateau, in one place reaching the opposite margin of the plateau surface. The magnitude of creep decreases steadily away from the head scarp (vi–x in Figs. 8c, d, and 9b). This zone of new deformation extends hundreds of metres down both lateral margins.

https://www.nat-hazards-earth-syst-sci.net/19/679/2019/nhess-19-679-2019-f09

Figure 9Displacement time series from the Pampahasi area (locations i to x in Fig. 8). (a) Movement of the west slope of the Río Irpavi valley within the area affected by the 2011 landslide (i–iii) and adjacent paleolandslide deposits (iv–v). The right axis is approximate, given that additional error is introduced from conversion from LOS to true motion. (b) Movement on Pampahasi Plateau behind the head scarp of the 2011 landslide (vi–x). Variability in the data, which is most pronounced in parts of the time series with low displacement rates, indicates phase noise resulting from factors unrelated to ground motion. Displacements between the end of the pre-failure stack and the start of the post-failure stack are not precisely constrained within the limits of the 2011 landslide (i–iii). For the Pampahasi slope adjacent to the 2011 failure and for the Pampahasi Plateau, displacements between the stacks are based on results of full-stack processing (see Sect. S4.8).

Download

All portions within the 2011 landslide area where coherence was maintained in both InSAR stacks and most areas within several hundred metres of the landslide experienced post-failure changes in slow ground motion (Fig. 10). The changes largely involved increased motion away from the satellite of between 2 and at least 5 cm a−1 LOS. Only a few small areas (totaling < 4 ha) in the northeast part of the toe, which were creeping at 4 cm a−1 LOS or more prior to February 2011, slowed following failure.

6 Conceptual failure model

Combining details of the 2011 failure event with patterns of long-term pre-failure and post-failure creep, we identify several types of instability that inform the geomechanical evolution of the Pampahasi slope.

6.1 Components of failure

The failing slope has many interconnected components that are largely within, but also extend locally beyond, the Pampahasi paleolandslide (Fig. 11a). Slope morphology before and after the 2011 landslide shows that long-term failure occurs predominantly by complex rotational–translational sliding (Fig. 11b). Differing displacement behaviour of several zones within the slope indicate separate, although likely interacting, kinematic components. The northern part of the Pampahasi slope toe consists of an apparently competent mass that creeped uniformly at rates of 2–5 cm a−1 LOS prior to late February 2011 (Fig. 8b) and moved more or less coherently during the 2011 landslide (Fig. 2b). This block-like mass accelerated slightly to 3–6 cm a−1 LOS following the 2011 failure and displayed a more spatially variable pattern of deformation (Fig. 8d) suggesting partial break-up during the failure event. The resulting removal of toe support – rapidly during the 2011 landslide but otherwise gradually – likely drove creep of the 12 ha zone of paleolandslide material (v) directly upslope (to the northwest). Although not part of a historic landslide event (Fig. 2), this zone was creeping slowly (∼0.5 cm a−1 LOS, Fig. 8a) prior to the 2011 failure and accelerated thereafter (> 1 cm a−1 LOS, Fig. 8c).

https://www.nat-hazards-earth-syst-sci.net/19/679/2019/nhess-19-679-2019-f10

Figure 10Change in line-of-sight slope creep of the Pampahasi area following the 2011 failure relative to pre-failure creep. Positive and negative values indicate, respectively, decreased and increased post-failure displacement rates away from the satellite. Green areas experienced no change in displacement and were largely stable both before and after the 2011 landslide. White areas lacked coherence phase in pre-failure scenes, post-failure scenes, or both. Roman numerals show locations of displacement histories in Fig. 9. Pre-failure and post-failure imagery spans, respectively, 8 September 2008 to 13 February 2011 and 9 March to 22 December 2011.

Download

https://www.nat-hazards-earth-syst-sci.net/19/679/2019/nhess-19-679-2019-f11

Figure 11Proposed geomechanical model for the Pampahasi landslide and adjacent terrain. (a) Components of failure, which differ in kinematics and activity history. (b) Generalized section through the failing slope illustrating the approximate depth and attitude of the failure zone. The topographic surface is based on NASA's Shuttle Radar Topographic Mission v3.0, which depicts the generalized pre-failure terrain. Dashed and solid lines in the source area (0 to ∼150 m x-axis distance) represent, respectively, pre-failure topography and post-failure topography based on terrain changes measured in the field using a handheld laser rangefinder. Terrain lower on the failed slope changed comparatively little during the landslide; the solid line beyond ∼150 m x-axis distance is thus generally representative of both the pre-failure and post-failure topography. (c) Schematic temporal displacement history generalized from InSAR measurements bracketing the 2011 Pampahasi landslide and eyewitness accounts. Due to the revisit frequency of RADARSAT-2 (24 days), the lack of acceleration leading up to the 2011 Pampahasi landslide is based on the lack of such effects reported by local residents. Velocity classes are from Cruden and Varnes (1996), spanning all classes except for “extremely rapid” (> 5 m s−1). The inset shows schematic stages of slope movement during first-time failure by fracture propagation (after Hermanns and Longva, 2012).

Download

A 1 ha zone near the midline of the landslide toe (ii) is one of the few locations of decreased post-failure activity (Fig. 10). Its high, localized activity prior to the 2011 landslide (5–8 cm a−1 LOS, Fig. 8b) suggests a shallow failure driven by incision at a location where Río Irpavi impinges on the toe of the slope (Fig. 2). This localized activity is probably similar to the east-bank landslide mapped by Anzoleaga et al. (1977) ∼400 m farther upstream, opposite the Río Jankopampa confluence (Fig. 2a), which was creeping throughout the entire period of radar monitoring (Fig. 8). Decreased post-failure activity of the small west-bank instability may reflect local modification of the Río Irpavi channel by the 2011 landslide and subsequent channel remediation efforts, or it may result from post-failure stress reorganization. Several similar linear areas of decreased post-failure creep along the toe of the 2011 landslide as far north as Río Jankopampa (Fig. 10) may represent similar conditions.

The southern part of the toe of the paleolandslide (iii in Fig. 8), which was stationary prior to the 2011 Pampahasi landslide (Fig. 8a), reactivated following the failure (Fig. 8c). Although the exact cause of the reactivation is uncertain, it is clearly a result of forces imposed on the toe of the paleolandslide by the 2011 event and may be related to consequent adjustment of stream erosion along the east margin of the deposit. Creep in this area was greatest directly upslope of the place where several buildings collapsed (Fig. 6f) in the final days (28 February and 1 March) of the landslide event, suggesting that delayed infrastructure damage there relates to the transition from the failure event to the new post-failure instability regime.

The height and steepness of the 2011 landslide head scarp indicate a deep-seated basal failure zone in the upper slope (Fig. 11b). Neither the paleolandslide nor the 2011 failure aligns with documented faults (Dobrovlny, 1962; Lavenu, 1978; Fig. 1b) that might indicate a structural control. A lack of displacement records for the head region following the 2011 landslide limits inferences about post-failure behaviour. However, pre-failure records, although incomplete, provide some insight into localized activity. Displacements within the head region prior to the landslide and the lack of movement in adjacent steep paleolandslide scarps suggest that pre-failure creep was part of the ongoing slope failure rather than localized instability resulting only from steep slopes. Localized small-scale slumping occurred lower in the head region along Río Chujilluncani (A in Fig. 8a, b) adjacent to, or potentially within, a region of pre-failure creep (∼1 cm a−1 LOS), supporting eyewitness accounts suggesting that the 2011 landslide initiated there.

The new displacements behind the head scarp and along the upper lateral margins of the 2011 landslide (Fig. 8c) gradually expanded over time into, respectively, the Pampahasi gravel (Fig. 9b) and paleolandslide deposits (iv and v in Fig. 9a). Movement in the latter suggests post-failure reactivation of paleolandslide material as a result of the 2011 landslide. Movements on the plateau are likely the result of dilation of the Pampahasi gravel in response to removal of material to a depth of up to 80 m along the 2011 head scarp (Fig. 6a).

6.2 Mechanisms of failure

Progressive failure is increasingly being recognized as a precursor of large, rapid landslides (e.g. Petley et al., 2002, 2005; Kilburn and Petley, 2003; Intrieri et al., 2018). In first-time brittle failures, creep reflects damage accumulation during the reduction of intact material strength, discontinuity strength, or both (Brideau and Roberts, 2015, and references therein) and consequent stress concentration (Cornelius and Scott, 1993). In such slopes, acceleration increases over time (Petley et al., 2002; Kilburn and Petley, 2003). Eventual development of a continuous failure surface (Petley et al., 2005) enables catastrophic release, commonly resulting in debris fragmentation, high run-out, and velocities exceeding 5 m s−1 (Hermanns and Longva, 2012). Sudden exceedance of resisting forces – for instance by seismic loading, increased pore pressure, or slope undercutting – may prematurely complete the failure surface, rapidly terminating the displacement trend (Hermanns and Longva, 2012). In slopes failing along a ductile shear zone or an existing plane of weakness, creep acceleration decreases over time, trending toward zero (Petley et al., 2002). Because displacement typically accelerates in each of the aforementioned failure types (Fig. 11c inset), it provides an opportunity for failure forecasting (Voight, 1989; Crosta and Agliardi, 2003; Eberhardt, 2012; Baron and Supper, 2013; Hermanns et al., 2013a; Federico et al., 2015; Moretto et al., 2017).

In contrast, over the 40-month period we monitored the Pampahasi slope (Fig. 11c), long-term steady-state creep was punctuated by temporary acceleration. Pre-2011 displacement rates are for the most part extremely slow (after Cruden and Varnes, 1996) but constant. No progressive acceleration was evident at the temporal resolution of RADARSAT-2 or from eyewitness accounts during the days leading up to failure. Displacement rates during the 4-day landslide event ranged from moderate to very rapid but did not reach the extremely rapid range (> 5 m s−1: Cruden and Varnes, 1996) typical of catastrophic failures. Creep rates following the 2011 landslide were also constant but were greater than those prior to failure, indicating a new state of instability following the landslide, which may reflect temporary post-failure adjustment. Constant creep rates suggest that the slope had, prior to failure and quickly afterward, reached steady-state behaviour, consistent with the presence of an existing failure surface (see Petley et al., 2002). Consequentially, prediction of large-scale failure events on the Pampahasi slope or at similar sites in La Paz should not rely exclusively on long-term creep time series.

6.3 Controlling factors

Geomechanical behaviour of the Pampahasi slope indicates failure nearer to that of sediment than to rock. Both the La Paz Formation and the Pampahasi gravel have compressive strengths well below 1 MPa, the threshold typically used to differentiate between rock and sediment (Brideau and Roberts, 2015). Parts of the main failure zone, as well as minor rupture zones within the landslide body, likely follow weak, fine-grained seams sheared during the Pampahasi paleolandslide or other past failures. Residual shear strength will be especially low in clay-rich zones. The lower half of the valley slope and the inferred failure zone beneath it dip at a similar angle (∼8, Fig. 11b) to the residual shear strength of the weakest (clay) zones within the failed sequence (Table S2), further suggesting that failure occurred along an existing failure surface.

Geotechnical properties of the slope, particularly its reduced strength resulting from extensive previous failure, enable its continual activity. However, external factors contribute to the long-term instability of this slope and many others in the La Paz basin. Ongoing fluvial downcutting and toe erosion by Río Irpavi help to maintain the slope's metastable state; spatiotemporal changes in fluvial erosion may have contributed to the differing pre-failure activity in the northern and southern parts of the landslide, as well as the increase in post-failure activity in both areas. This driver is a response to the 4 km drop in base level since the Early Pleistocene – from the time when the Altiplano plateau was internally draining to the present drainage to the Atlantic Ocean (Roberts et al., 2018) – and operates along other trunk streams in the basin.

Fewer than half of the largest historic landslides in the La Paz area follow intense precipitation events (Table S1), suggesting that enhanced pore-water pressure may play a role in initiating only some large failures. The sudden onset of the 2011 Pampahasi landslide (Fig. 11c) is typical of triggered landslides (Fig. 11c inset; Hermanns and Longva, 2012). Heavy precipitation on the previous day – the 10th wettest on record – is the only apparent trigger for the failure. Slope wetting was likely enhanced by effluent from homes and concentrated overland flow from impervious surfaces. We thus postulate that rainwater infiltrating the slope and plateau above it on 25 February reached the basal failure zone by the following evening; a consequent rapid pore-pressure increase and shear-strength reduction appear to have triggered the landslide. Rupture of the water line crossing the slope added much more water and likely contributed to the mobility of the 2011 landslide directly downslope, where displacement rates were greatest.

6.4 Sources of uncertainty

Possible acceleration in the days or weeks prior to the 2011 landslide cannot be evaluated with the temporal resolution of the displacement histories used here. RADARSAT-2's 24-day revisit period provides low sampling frequency compared to some other radar satellites (Moretto et al., 2017) and especially to in situ sensors or ground-based remote sensors (e.g. Kalaugher et al., 2000; Crosta and Agliardi, 2003; Gischig et al., 2011; Eberhardt, 2012; Federico et al., 2015; Confuorto et al., 2017; Carlà et al., 2018). Additionally, the temporal filter used in processing (Sect. S4.6) attenuates records of non-linear displacement trends. The resulting temporal details reported here are sufficient to characterize gradual acceleration over periods of months (Fig. 9), as has been documented for pre-failure creep using InSAR elsewhere (Carlà et al., 2018; Intrieri et al., 2018), or more abrupt acceleration over just a few radar acquisitions but not shorter-term changes. Sporadic displacement activity that may signal impending acceleration (see Kalaugher et al., 2000) could similarly have gone undocumented.

Ground motion represented in displacement maps is independent of the structural behaviour of the built environment. Isolated building instability is likely in light of some local construction practices in La Paz but will be extremely localized and thus is removed during spatial interpolation of the maps. Phase change due to thermal expansion will be minimal given limited seasonal temperature differences in the study area. Due to their cyclic nature, any such phase component will not influence long-term displacement trends.

The HDS-InSAR processing chain is complex and includes many non-linear steps, which greatly complicates development of an accurate error-approximation model. Both the pre-failure and post-failure stacks have good baseline diversity, allowing relative errors between them to be approximated by first-order estimates from the square root of the number of scenes (32 vs. 12). In the absence of a rigorous model, we assume that error for the thinner stack is conservatively twice as large as that of the thicker stack. We approximate the errors as 3 and 6 mm a−1, respectively, for the pre-failure and post-failure stacks. Because the average displacement rates across much of the Pampahasi area exceed the error estimates, the displacement patterns documented here are deemed to be reliable.

Due to the structure of the HDS-InSAR processing chain, differing environmental conditions between the two stacks, namely topography and moisture (see Wasowski and Bovenga, 2014, 2015), have minimal effects. The reference digital terrain model (DTM) is used only for an initial topographic correction; stack processing solves for height error relative to the reference DTM and provides a new elevation solution for each of the two stacks (Sect. S4.6), which improves terrain representation. Topographic correction of the post-failure stack thus accounts for landslide-induced terrain changes, which were greatest in the source area (Fig. 11b). Temporal soil moisture variability is unlikely to affect phase by more than 100 (Rabus et al., 2010), which equates to approximately one-third of an interferometric fringe or 0.9 cm for the sensor used here. Comparison of precipitation records in the 30 months before and 10 months after the 2011 Pampahasi landslide indicates that long-term precipitation amounts during the pre-failure and post-failure stacks were comparable (Fig. S2; Table S4). Spatial moisture gradients are a more substantial error source (Rabus et al., 2010), but major differences present in a single scene are removed during stack processing.

The assumption of slope-parallel motion is likely incorrect in some parts of the Pampahasi area. The complex nature of the 2011 failure, particularly the variability in the displacement directions (Fig. 2b), suggests that separate masses may have moved obliquely relative to one another during pre-failure and post-failure creep. Localized compression or extension within the landslide would additionally impart vertical components of movement. Additional error is thus introduced in the conversion from measured one-dimensional (LOS) displacement to approximate three-dimensional (true) motion. Provided that the sense of movement at any given location on the slope has not substantially changed following the 2011 failure, LOS motion vectors record similarly scaled changes in the true rate of motion and thus reliably represent the degree of change in slope activity. Additionally, the close alignment between the radar LOS and the fall line of the slope ensures that the InSAR-measured displacement record is particularly sensitive to expected gravity-driven slope motion between Pampahasi Plateau and Río Irpavi.

Because RADARSAT-2 acquisitions over the site ended 10 months after the failure, the current activity of the slope is poorly constrained and it is possible that increased post-failure slope creep was short-lived. Additional InSAR acquisitions planned for this and other sites in La Paz, including both ascending and descending orbits, will provide additional insight into more recent temporal and spatial complexity of slope activity.

7 Advances in InSAR monitoring of landslides

7.1 Comparison with other studies

Our detailed documentation of slope activity preceding and following a large landslide using spaceborne InSAR is, to our knowledge, unique. Numerous studies show evidence of long-term, generally sustained creep of slow-moving landslides (e.g. Colesanti and Wasowski, 2006; Lauknes et al., 2010; Hermanns et al., 2013b; Necsoiu et al., 2014; Schlögel et al., 2016) and in many cases document temporal variability in displacements. A much smaller number of studies identify short-term deformation between successive scenes (Moro et al., 2011) or creep leading up to catastrophic failure (Carlà et al., 2018; Intrieri et al., 2018; Handwerger et al., 2019). However, we are aware of only one other spaceborne InSAR investigation that compares pre-failure and post-failure activity. Confuorto et al. (2017) report displacement rates during a 1-year period starting ∼20 months following the 2012 Via Piave landslide in southern Italy that were marginally faster than rates during a ca. 2-year period ending ∼16 months before failure (see Sect. S5 for comparison with the dataset presented here). In that investigation, the limited number of point reflectors, their concentration in a narrow band of the head scarp, and the 3-year data gap bracketing the landslide limit insight into potential post-failure changes in the slope.

7.2 Expanding utility

Advanced InSAR approaches such as HDS-InSAR increase ground-target density and discriminate spatially uncorrelated deformation, improving the characterization of landslides with spaceborne InSAR. Even though applied here to radar data of moderate spatial and temporal resolution, HDS-InSAR has provided detailed characterization of a large, dynamic urban slope. Applying the technique to datasets with higher spatial resolution will allow even more detailed characterization of landslide creep, and sensor constellations will help increase the temporal resolution of displacement records. The HDS technique also provides insight into the precise locations of landslide margins and is thus particularly useful in dense urban settings where development and remedial works remove evidence that is useful in landslide investigations.

Furthermore, such techniques enable more comprehensive detection and characterization of landslides of different depth and size. Preferential detection of deep landslides by InSAR (Wasowski and Bovenga, 2014, 2015) reflects the typically high spatial regularity of their displacements. Shallow instability, especially in areas of variable microtopography, generates spatially irregular ground motion that is more difficult to detect. HDS-InSAR's optimization for locally variable ground motion particularly improves the characterization of shallow landslides and thus reduces biasing toward deep-seated instability.

8 Conclusions

The 2011 Pampahasi landslide is only one event in the long-term, and likely ongoing, activity of a metastable slope. High-spatial-density InSAR results demonstrate that several parts of the slope were active over the years prior to the landslide and that nearly all regions with coherent ground-motion records within and adjacent to it experienced post-failure increases in the extent and rate of creep. This change in ground deformation counters any expectation that a complex landslide might stabilize, at least temporarily, following a discrete failure and highlights that such an assumption, at least for short-term stability, is imprudent for multigenerational failures. Slope dynamics documented here support observations and theory from the scientific literature that spatiotemporal clustering of landslides are responses to stress redistribution, which is likely to exceed stress changes due to background erosion.

InSAR-derived deformation patterns provide new insight into failure mechanisms of the Pampahasi slope. The complex failure comprises multiple unstable parts, including many areas of shallow instability and at least one large creeping block. Accelerating creep, a typical feature of first-time and reactivated failures, was not apparent leading up to the disastrous landslide on 26 February 2011. Any increase in displacement during the 13-day window following the last pre-failure radar scene (13 February) apparently was not perceived by residents on the slope. The absence of detectable acceleration may instead reflect the achievement of steady-state creep along an existing surface, which is supported by the presence of previously failed sediments and the similarity of the gradient of lower slopes to the residual strength of the weakest geologic units. Increased pore pressures resulting from particularly heavy rainfall triggered a large failure that temporarily perturbed the dynamic equilibrium of the slope. The failure reorganized stresses in the slope as well as beneath the plateau above it, leading to a new, although possibly temporary, state of dynamic equilibrium affecting even more of the slope.

Improved understanding of instability of the Pampahasi slope is instructive in evaluating and reducing risk from landslides in La Paz, especially given the large extent of paleolandslides in the city, their close association with recent landslides, and incomplete knowledge of slope activity. The Pampahasi slope is similar to numerous other slopes formed on large ancient landslides in the city, including slopes where many of the largest historic landslides have occurred (Fig. 1b). These slopes are underlain by generally weak, fine-grained sediments of the La Paz Formation. Such slopes, particularly those descending from the west margin of Pampahasi Plateau, are key sites for future stability and risk assessments. Additionally, stress redistribution suggested by enhanced creep of the Pampahasi paleolandslide between the margins of the 2009 and 2011 failures and possibly in adjacent parts of the Villa Salomé paleolandslide may increase their susceptibility to large-scale rapid failure in the future. Evaluating the possibility of future large-scale reactivations of complex landslides should not, however, rely solely on creep acceleration, especially given the lack of such a signal leading up to the 2011 landslide. Furthermore, although fewer than half of the historic failures exceeding 1 Mm3 happened during the rainy season, coincidence of the 2011 reactivation of the Pampahasi paleolandslide with particularly wet conditions indicates that consideration of high-rainfall scenarios is advisable. In light of the evidence presented here of enhanced post-failure activity, creeping slopes in La Paz are generally not appropriate for new settlements or resettlement of residents displaced by disasters.

Optimizing InSAR characterization of landslides requires processing methodologies that afford high ground-target density and discriminate spatially uncorrelated deformation. The improved understanding of the activity of the Pampahasi area reflects the detailed spatial characterization made possible by the HDS-InSAR technique. Applying advanced InSAR processing techniques such as this to modern and next-generation radar datasets that, due to continued technological improvement, have higher ground resolution and revisit frequency will improve understanding of complex unstable slopes.

Data availability

The data are not publicly available as HDS-InSAR is part of a commercial product line. The processed data contain intellectual property developed through research and development at MacDonald, Dettwiler and Associates (MDA).

Supplement

The supplement related to this article is available online at: https://doi.org/10.5194/nhess-19-679-2019-supplement.

Author contributions

NJR, JJC, RLH, and BTR conceived the investigation. NJR processed radar datasets with guidance from BTR. All authors participated in fieldwork. MAG and EM arranged field logistics and provided crucial details on the history of development and instability. NJR, JJC, BTR, and RLH interpreted observations and results. NJR prepared the paper with contributions from all co-authors.

Competing interests

The authors declare that they have no conflict of interest.

Acknowledgements

Nicholas J. Roberts used facilities and software at MacDonald, Dettwiler and Associates (MDA) courtesy of Harold Zwick and Christian Nadeau and was assisted in house by Jayanti Sharma and Jayson Eppler. He was partially funded by the Natural Science and Engineering Research Council (NSERC) of Canada (postgraduate scholarship) and the American Society for Photogrammetry and Remote Sensing (ASPRS; Robert C. Colwell Fellowship). MDA Geospatial Services provided RADARSAT-2 imagery. Field investigation was supported by NSERC (Discovery Grant 24595 to John J. Clague), APSRS (Ta Liang Memorial Award to Nicholas J. Roberts), and the International Centre for Geohazards (ICG contribution no. 377 to Reginald L. Hermanns). Victor Ramírez and Eddy Baldellón (municipality of La Paz) relayed eyewitness accounts of the 2011 Pampahasi landslide reported by residents of the affected area. Detailed comments and suggestions from the anonymous reviewers significantly improved the final paper.

Edited by: Thomas Glade
Reviewed by: two anonymous referees

References

Aguilar, O. Q.: Detección de cambios de uso del suelo en la zona de Callapa afectada por el megadeslizamiento (La Paz – Bolivia), in: Proceedings of the XIV Encuentro de Geógrafos de América Latina, Lima, Peru, 8–11 April 2013, Encuentro de Geógrafos de América Latina, Lima, Peru, 20 pp., 2013. 

Ahlfeld, F.: Reseña geológica de la Cuenca de La Paz – Parte 1, Instituto Boliviano de Ingeniería de Minas y Geología, 2, 11–15, 1945a. 

Ahlfeld, F.: Reseña geológica de la Cuenca de La Paz – Parte 2, Instituto Boliviano de Ingeniería de Minas y Geología, 2, 11–15, 1945b. 

Anzoleaga, O., Malatrait, A. M., Michalski, E. R., Bles, J. L., Otazo, N., and Walgenwitz, F.: Características geomtécnicas de la cuenca de La Paz y alrededores, Plan de desarrollo Urbano para la ciudad de La Paz, Honorable Alcaldía Municipal de La Paz, Bolvia, Informe Geológico no. 15, 107 pp., 1977. 

Baron, I. and Supper, R.: Application and reliability of techniques for landslide site investigation, monitoring and early warning – outcomes from a questionnaire study, Nat. Hazards Earth Syst. Sci., 13, 3157–3168, https://doi.org/10.5194/nhess-13-3157-2013, 2013. 

Berardino, P., Fornaro, G., Lanari, R., and Sansosti, E.: A new algorithm for surface deformation monitoring based on small baseline differential SAR interferograms, IEEE T. Geosci. Remote, 40, 2375–2383, https://doi.org/10.1109/TGRS.2002.803792, 2002. 

Bles, J. L., Alvarez, A., Anzoleaga, O., Ballivián, O., Bustillos, O., Hochstatter, H., Malatrait, A. M., and Otazo, N.: Características litoestratigráficas de la cuenca de La Paz y alrededores, Plan de desarrollo Urbano para la ciudad de La Paz, Honorable Alcaldía Municipal de La Paz, Bolivia, Informe Geológico no. 5, 35 pp., 1977. 

Brideau, M. A. and Roberts, N. J.: Mass movements in bedrock, in: Landslides Hazards, Risks and Disasters, edited by: Davies, T. and Shroder, J. F., Academic Press, Amsterdam, Netherlands, 43–90, https://doi.org/10.1016/B978-0-12-396452-6.00003-3, 2015. 

Carlà, T., Farina, P., Intrieri, E., Ketizmen, H., and Casagli, N.: Integration of ground-based radar and satellite InSAR data for the analysis of an unexpected slope failure in an open-pit mine, Eng. Geol., 235, 39–52, https://doi.org/10.1016/j.enggeo.2018.01.021, 2018. 

Chigira, M. and Kiho, K.: Deep-seated rockslide–avalanches preceded by mass rock creep of sedimentary rocks in the Akaishi Mountains, Central Japan, Eng. Geol., 38, 221–230, https://doi.org/10.1016/0013-7952(94)90039-6, 1994. 

Colesanti, C. and Wasowski, J.: Investigating landslides with space-borne Synthetic Aperture Radar (SAR) interferometry, Eng. Geol., 88, 173–199, https://doi.org/10.1016/j.enggeo.2006.09.013, 2006. 

Confuorto, P., Di Martire, D., Centolanza, G., Iglesias, R., Mallorqui, J. J., Novellino, A., Plank, S., Ramondini, M., Thuro, K., and Calcaterra, D.: Post-failure evolution analysis of a rainfall-triggered landslide by multi-temporal interferometry SAR approaches integrated with geotechnical analysis, Remote Sens. Environ., 188, 51–72, https://doi.org/10.1016/j.rse.2016.11.002, 2017. 

Cornelius, R. R. and Scott, P. A.: A materials failure relation of accelerating creep as empirical description of damage accumulation, Rock Mech. Rock Eng., 26, 233–252, https://doi.org/10.1007/BF01040117, 1993. 

Crosta, G. B. and Agliardi, F.: Failure forecast for large rock slides by surface displacement measurements, Can. Geotech. J., 40, 176–191, https://doi.org/10.1139/t02-085, 2003. 

Crosta, G. B., Hermanns, R. L., Dehls, J., Lari, S., and Sepulveda, S.: Rock avalanches clusters along the northern Chile coastal scarp, Geomorphology, 289, 27–43, https://doi.org/10.1016/j.geomorph.2016.11.024, 2017. 

Cruden, D. M. and Varnes, D. J.: Landslide types and processes, in: Landslides Investigation and Mitigation: Spec. Rep. 247, edited by: Turner, A. K. and Schuster, L. R., Transportation Res. Board, Washington, D.C., 36–75, 1996. 

Dobrovolny, E.: Preliminary report on the geology of the upper La Paz Valley, U.S. Geol. Surv., Open File Rep., 55–43, 32 pp., 1955. 

Dobrovolny, E.: Geología del valle de La Paz, Departamento Nacional de Geología del Ministerio de Minas y Petroleo, Servicio Geológico de Bolivia Boletín 3, 152 pp., 1962. 

Dobrovolny, E.: A postglacial mudflow of large volume in the La Paz valley, Bolivia, U.S. Geol. Surv. Prof. Pap., 600-C, 130–134, 1968. 

Eberhardt, E.: Landslide monitoring: The role of investigative monitoring to improve understanding and early warning of failure, in: Landslides: Types, Mechanisms and Modeling, edited by: Clague, J. J. and Stead, D., Cambridge University Press, Cambridge, UK, 222–234, https://doi.org/10.1017/CBO9780511740367.020, 2012. 

Eppler, J. and Rabus, B.: Monitoring urban infrastructure with an adaptive multilooking InSAR technique, in: Proceedings of the Fringe Workshop, Frascati, Italy, 19–26 September 2011, European Space Agency Communications SP-697, Unpagninated CD ROM, 5 pp., 2011. 

Federico, A., Popescu, M., and Murianni, A.: Temporal prediction of landslide occurrence: A possibility or a challenge?, Italian J. Eng. Geol. Environ., 1, 41–60, https://doi.org/10.4408/IJEGE.2015-01.O-04, 2015. 

Ferretti, A., Prati, C., and Rocca, F.: Permanent scatterers in SAR interferometry, IEEE T. Geosci. Remote, 39, 8–20, https://doi.org/10.1109/36.898661, 2001. 

Ferretti, A., Fumagalli, A., Novali, F., Prati, C., Rocca, F., and Rucci, A.: A new algorithm for processing interferometric data-stacks – SqueeSAR, IEEE T. Geosci. Remote, 49, 3460–3470, https://doi.org/10.1109/TGRS.2011.2124465, 2011. 

Gischig, V., Amann, F., Moore, J. R., Loew, S., Eisenbreiss, H., and Stempfhuber, W.: Composite rock slope kinematics at the current Randa instability, Switzerland, based on remote sensing and numerical modeling, Eng. Geol., 118, 37–53, https://doi.org/10.1016/j.enggeo.2010.11.006, 2011. 

Handwerger, A. L., Huang, M.-H., Fielding, E. J., Booth, A. M., and Bürgmann, R.: A shift from drought to extreme rainfall drives a stable landslide to catastrophic failure, Sci. Rep., 9, 1569, https://doi.org/10.1038/s41598-018-38300-0 1, 2019. 

Hermanns, R. L. and Longva, O.: Rapid rock-slope failures, in: Landslides: Types, Mechanisms and Modeling, edited by: Clague, J. J. and Stead, D., Cambridge Univ. Press, Cambridge, UK, 59–70, https://doi.org/10.1017/CBO9780511740367.007, 2012. 

Hermanns, R. L., Blikra, L. H., Naumann, M., Nilsen, B., Panthi, K. K., Stromeyer, D., and Longva, O.: Examples of multiple rock-slope collapses from Köfels (Ötz valley, Austria) and western Norway, Eng. Geol., 83, 94–108, https://doi.org/10.1016/j.enggeo.2005.06.026, 2006. 

Hermanns, R. L., Dehls, J. F., Guzmán, M.-A., Roberts, N. J., Clague, J. J., Cazas Saavedra, A., and Quenta, G.: Relation of recent megalandslides to prehistoric events in the city of La Paz, Bolivia, in: Proceedings of the 2nd North American Symposium on Landslides, 3–8 June 2012, Banff, Alberta, Volume 1, Taylor & Francis Group, London, UK, 341–347, 2012. 

Hermanns, R. L., Oppikofer, T., Anda, E., Blikra, L. H., Böhme, M., Bunkholt, H., Crosta, G. B., Dahle, H., Devoli, G., Fischer, L., Jaboyedoff, M., Loew, S., Sætre, S., and Yugsi Molina, F.: Hazard and risk classification system for large unstable rock slopes in Norway, in: International Conference on Vajont – 1963–2013, edited by: Genevois, R. and Prestininzi, A., Italian J. Eng. Geol. Environ., Book Ser., 6, 245–254, https://doi.org/10.4408/IJEGE.2013-06.B-22, 2013a. 

Hermanns, R. L., Oppikofer, T., Dahle, H., Eiken, T., Ivy-Ochs, S., and Blikra, L. H.: Understanding long-term slope deformation for stability assessment of rock slopes: The case of the Oppstadhornet rockslide, Norway, in: International Conference on Vajont – 1963–2013, edited by: Genevois, R. and Prestininzi, A., Italian J. Eng. Geol. Environ., Book Ser., 6, 255–264, https://doi.org/10.4408/IJEGE.2013-06.B-23, 2013b. 

Hilger, P., Hermanns, R. L., Gosse, J. C., Jacobs, B., Etzelmüller, B., and Krautblatter, M.: Multiple rock-slope failures from Mannen in Romsdal Valley, western Norway, revealed from Quaternary geological mapping and 10Be exposure dating, Holocene, 28, 1841–1854, https://doi.org/10.1177/0959683618798165, 2018. 

Intrieri, E., Raspini, F., Fumagalli, A., Lu, P., Del Conte, S., Farina, P., Allievi, J., Ferretti, A., and Casagli, N.: The Maoxian landslide as seen from space: Detecting precursors of failure with Sentinel-1 data, Landslides, 15, 123–133, https://doi.org/10.1007/s10346-017-0915-7, 2018. 

Kalaugher, P. G., Hodgson, R. L. P., and Granger, P.: Pre-failure strains as precursors of sliding in a coastal mudslide, Q. J. Eng. Geol. Hydroge., 33, 325–334, https://doi.org/10.1144/qjegh.33.4.325, 2000. 

Kilburn, C. R. J. and Petley, D. N.: Forecasting giant, catastrophic slope collapse; lessons from Vajont, Geomorphology, 54, 21–32, https://doi.org/10.1016/S0169-555X(03)00052-7, 2003. 

Lauknes, T. R., Piyush Shanker, A., Dehls, J. F., Zebker, H. A., Henderson, I. H. C., and Larsen, Y.: Detailed rockslide mapping in northern Norway with small baseline and persistent scatterer interferometric SAR time series methods, Remote Sens. Environ., 114, 2097–2109, https://doi.org/10.1016/j.rse.2010.04.015, 2010. 

Lavenu, A.: Néotectonique des sediments plio-quaternaires de l'Altiplano bolivien (region de La Paz, Ayo Ayo, Umala), Cah. ORSTOM Géol., 10, 115–126, 1978. 

Lavenu, A., Bonhomme, M. G., Vatin-Perignon, N., and de Pachtere, P.: Neogene magmatism in the Bolivian Andes between 16 S and 18 S: Stratigraphy and K/Ar geochronology, J. S. Am. Earth Sci., 2, 35–47, https://doi.org/10.1016/0895-9811(89)90025-4, 1989. 

Moretto, S., Bozzano, F., Esposito, C., Mazzanti, P., and Rocca, A.: Assessment of landslide pre-failure monitoring and forecasting using satellite SAR interferometry, Geosciences, 7, 7020036, https://doi.org/10.3390/geosciences7020036, 2017. 

Moro, M., Chini, M., Saroli, M., Atzori, S., Stramondo, S., and Salvi, S.: Analysis of large, seismically induced, gravitational deformations imaged by high-resolution COSMO-SkyMed synthetic aperture radar, Geology, 39, 527–530, https://doi.org/10.1130/G31748.1, 2011. 

Necsoiu, M., McGinnis, R. N., and Hooper, D.M.: New insights on the Salmon Falls Creek Canyon landslide complex based on geomorphological analysis and multitemporal satellite InSAR techniques, Landslides, 11, 1141–1153, https://doi.org/10.1007/s10346-014-0523-8, 2014. 

O'Hare, G. and Rivas, S.: The landslide hazard and human vulnerability in La Paz City, Bolivia, Geogr. J., 171, 239–258, https://doi.org/10.1111/j.1475-4959.2005.00163.x, 2005. 

Petley, D. N., Bulmer, M. K., and Murphy, W.: Patterns of movement in rotational and translational landslides, Geology, 30, 719–722, https://doi.org/10.1130/0091-7613(2002)030<0719:POMIRA>2.0.CO;2, 2002. 

Petley, D. N., Higuchi, T., Petley, D. J., Bulmer, M. H., and Carey, J.: Development of progressive landslide failure in cohesive materials, Geology, 33, 201–204, https://doi.org/10.1130/G21147.1, 2005. 

Quenta, G. and Calle, A.: Diagnostico de estabilidad de la zona Valle de las Flores, [Proyecto Multinacional Andino: Geociencias para las Comunidades Andinas], Servicio Nacional de Geología y Tecnico de Minas, La Paz, Bolivia, 12 pp., 2005. 

Quenta, G., Galarza, I., Teran, N., Hermanns, R. L., Cazaz, A., and García, H.: Deslizamiento traslacional y represamiento en el valle de Allpacoma, ciudad de La Paz, Bolivia, in: Movimientos en Masa en la Región Andina – Una Guía para la Evaluación de Amenazas, Proyecto Multinacional Andino: Geociencias para las Comunidades Andinas, Servicio Nacional de Geología y Minería, Santiago, Chile, Publicación Geológica Multinacional, no. 4, 230–234, 2007. 

Quenta, G., Hermanns, R., Galarza, I., Teran, N., Cazas, A., Caballero, T., and García, H.: Deslizamiento y represamiento en el valle de Allpacoma, in: Experiencias Andinas en Mitigación de Riesgos Geológicos, Proyecto Multinacional Andino: Geociencias para las Comunidades Andinas, Servicio Nacional de Geología y Minería, Santiago, Chile, Publicación Geológica Multinacional, no. 6, 295–2833, 2008.  

Rabus, B., When, H., and Nolan, M.: The importance of soil moisture and soil structure for InSAR phase and backscatter, as determined by FDTD modeling, IEEE T. Geosci. Remote, 48, 2421–2429, https://doi.org/10.1109/TGRS.2009.2039353, 2010. 

Rabus, B., Eppler, J., Sharma J., and Busler, J.: Tunnel monitoring with an advanced InSAR technique, in: Proceedings of RADAR Sensor Technology XVI, 23–25 April 2012, Baltimore, USA, SPIE Press, Bellingham, WA, 10 pp., https://doi.org/10.1117/12.918644, 2012. 

Roberts, N. J.: Late Cenozoic geology of La Paz, Bolivia, and its relation to landslide activity, PhD thesis, Department of Earth Sciences, Simon Fraser University, Canada, 248 pp., 2016. 

Roberts, N. J., Barendregt, R. W., and Clague, J. J.: Multiple tropical Andean glaciations during a period of late Pliocene warmth, Sci. Rep., 7, 41878, https://doi.org/10.1038/srep41878, 2017. 

Roberts, N. J., Barendregt, R. W., and Clague, J. J.: Pliocene and Pleistocene chronostratigraphy of continental sediments underlying the Altiplano at La Paz, Bolivia, Quarternary Sci. Rev., 189, 105–126, https://doi.org/10.1016/j.quascirev.2018.03.008, 2018. 

Scanvic, J.-Y. and Girault, F.: Imagerie SPOT-1 et inventaire des mouvements de terrain: l'exemple de La Paz (Bolivia), Photo Interpretation, 28, 1–20, 1989. 

Schlögel, R., Malet, J.-P., Doubre, C., and Lebourg, T.: Structural control on the kinematics of the deep-seated La Clapière landslide revealed by L-band InSAR observations, Landslides, 13, 1005–1018, https://doi.org/10.1007/s10346-015-0623-0, 2016. 

Servant, M., Sempere, T., Argollo, J., Bernat, M., Ferand, G., and Lo Bello, P.: Morphogenèse et soulèvement de la Cordillère Orientale des Andes de Bolivie au Cénozoïque, CR Acad. Sci., 309, 417–422, 1989. 

Voight, B.: A relationship to describe rate-dependant material failure, Science, 243, 200–203, https://doi.org/10.1126/science.243.4888.200, 1989. 

Wasowski, J. and Bovenga, F.: Investigating landslides and unstable slopes with satellite Multi Temporal Interferometry: Current issues and future perspectives, Eng. Geol., 174, 103–138, https://doi.org/10.1016/j.enggeo.2014.03.003, 2014. 

Wasowski, J. and Bovenga, F.: Remote sensing of landslide motion with emphasis on satellite multitemporal interferometry applications: an overview, in: Landslides Hazards, Risks and Disasters, edited by: Davies, T. and Shroder, J. F., Academic Press, Amsterdam, Netherlands, 345–403, https://doi.org/10.1016/B978-0-12-396452-6.00011-2, 2015. 

1

The 2011 failure is also locally called the Callapa landslide after one of the neighbourhoods it destroyed.

2

This stream is sometimes also referred to as Arroyo Pampahasi.

Download
Short summary
La Paz, Bolivia, experiences frequent damaging landslides. We quantify creep before and after the city’s largest modern landslide using spaceborne InSAR. Creep of ancient landslide deposits increased in rate and extent following failure and extended into adjacent intact materials. Accelerated steady-state creep reflects complex post-failure stress redistribution. Landslide risk in La Paz, which is underlain by many large ancient landslides, may be even greater than previously thought.
Altmetrics
Final-revised paper
Preprint