Adapting sudden landslide identification product (SLIP) and detecting real-time increased precipitation (DRIP) algorithms to map rainfall-triggered landslides in Western Cameroon highlands (Central-Africa)

NASA’s developers recently proposed the Sudden Landslide Identification Product (SLIP) and Detecting Real-Time Increased Precipitation (DRIP) algorithms. This double method uses Landsat 8 satellite images and daily rainfall data for a real-time mapping of this geohazard. This study adapts the processing to face the issues of data quality and unavailability/gaps for the mapping of the recent landslide events in west-Cameroon’s highlands. The SLIP algorithm is adapted, by integrating the inverse Normalized Difference Vegetation Index (NDVI) to assess the soil bareness, the Modified Normalized Multi-Band Drought Index (MNMDI) combined with the hydrothermal index to assess soil moisture, and the slope inclination to map the recent landslide. Further, the DRIP algorithm uses the mean daily rainfall to assess the thresholds corresponding to the recent landslide events. Their probability density function (PDF) curves are superimposed and their intersections are used to propose sets of dichotomous variables before (1948–2018) and after the 28 October 2019 landslide event. In addition, a survival analysis is performed to correlate landslide occurrence to rainfall, with the first known event in Cameroon as starting point, and using the Cox model. From the SLIP model, the Landslide Hazard Zonation (LHZ) map gives an overall accuracy of 96%. Further, the DRIP model states that 6/9 ranges of probability are rainfall-triggered landslides at 99.99%, between June and October, while 3/9 ranges show only 4.88% of risk for the same interval. Finally, the survival probability for a known site is up to 0.68 for the best value and between 0.38 and 0.1 for the lowest value through time. The proposed approach is an alternative based on data (un)availability, completed by the site’s lifetime analysis for a more flexibility in observation and prediction thresholding.


Background
Landslides are natural events (Varnes, 1958(Varnes, , 1978(Varnes, & 1984Brusden, 1984;Crozier, 1986;Hutchinson, 1988;Cruden, 1991;Cruden and Varnes, 1996;Courture, 2011;UNISDR, 1 2017;USGS * , 2004). However, they may turn into serious geohazards responsible for casualties and economical losses worldwide (Petley, 2012). These include loss of lives and damage to human settlements and natural structures, which present a significant constraint for the development of the zones affected. According to the World Health Organization, landslides affected an estimated 4.8 million people and caused more than 18,000 deaths, between 1998 and 2017. 2 It is admitted that at least 90% of losses related can be avoided if the problem is recognized before the landslide occurrence (Brabb, 1993). Therefore, the mitigation measures require to identify existing landslides, and/or to predict of future events and endangered zones. One main issue is that landslide inventories suffer from underreporting at both regional and global scales (UNISDR, 2017). Even in developed countries, the database of landslide events is usually far from complete. Significant gaps in available information additionally contribute to the shortcomings of the inventories due to the lack of routine global monitoring or cataloging systems, such as those available for hurricanes and earthquakes (Kirschbaum et al. 2009). However, mapping landslide deformation and activity is fundamental for the assessment and reduction of hazards and risks related (Zhao and Lu, 2018). To the measures and sketches that result from fieldwork (Yang et al., 2012(Yang et al., & 2015, Earth Observatory (OE) brings the above view that is complementary to assess or predict the hazard.
Remote sensing data and the geospatial sciences are very powerful tools to study the prevailing causal factors and achieve that goal (Tofani et al., 2013). Their integration leads to a standard tool known as landslide susceptibility mapping used around the world by various researchers (Guzzetti et al., 1999;Van Westen et al., 1999), which helps mapping the areas affected or model the likelihood of future landslides based solely on the intrinsic properties of a site. The susceptibility of a given area to landslides can be determined and depicted using hazard zonation (Lin et al., 2017;Dahal and Dahal, 2017). Recent advances identify two sets of methods for landslide hazard zonation (LHZ), such as heuristic (knowledge-based) and data driven (statistical) approaches (Pardeshi et al., 2013;Roy and Saha, 2019).
Besides, a novelty approach has recently emerged, ambitioning to fill the methodological voids of landslide monitoring. In a double set of processing, based on the availability of hourly rainfall data and the more recent free satellite images, instantaneous assessment and prediction of landslide events became possible for spatially extended regions thanks to the analytical power and the flexibility of Google Earth Engine platform and tools. The pioneering model was developed in 2019 by a team of NASA's developers (Fayne et al., 2019) and is detailed in the methodology section. This method is known as, the Sudden Landslide Identification Product (SLIP), combined with the Detecting Real-Time Increased Precipitation (DRIP), to simultaneously map the affected surface, as well as identifies intensity and timing of rainfall that had triggered the event.
With the advantages of this latter, i.e., the opensource data and the relative flexibility of the method that relies on indices calculation and rainfall sequences, the present study has envisaged to map the rainfall-triggered landslides events in west-Cameroon highlands, based on recent deadly events. Nonetheless, there are several local obstacles to its full implementation, among which three majors: i) first of all, free Landsat 8 images covering the area usually have at least 50% cloud cover in permanence, and those less cloudy are unequally distributed in time, thus making impossible an efficient monitoring every 16 days; ii) the poor spatial and temporal data availability, due to the insufficient number of stations only present in few main cities, and the lack of hourly records caused by obsolete materials; iii) at last, and non the least issue, there is no a detailed and accurate landslide inventory in Cameroon, in terms of date, proper location with geographical coordinates, and statistics on social and economic losses. These lacks and gaps clearly impede the full implementation of the SLIP/DRIP model. Therefore, this paper proposes a derivative SLIP/DRIP procedure, adapted to countries and world's areas facing the same difficulties above-listed. This alternative fits in this study, Cameroon's local scientific and economic conditions, by finding way-outs for efficient results, adjusting the processing as possible to a context of missing data, and introducing a complementary process of the survival analysis of (potentially) affected sites.

Study location
This study was carried out at a regional scale, on a subset of Cameroon's western highlands (Fig. 1a-c) covering 3930.35 km 2 . The area belongs to the Cameroon Volcanic Line (CVL), one of the several segments of the African Plate, oriented NE between 9 and 11°of latitude and SW 5-7°East longitude and forming a horst (Elsheikh et al., 2014;. It is an area of transition between the Cameroon's rainy and dry areas. Its morphology is complex and consists of high plateaus, volcanic massifs as well as plains or collapse basins corresponding to the graben (Balla et al., 2013). The three main morphological components are the Bamenda Plateau (2200 m), the Bamiléké Plateau (1400 to 1600 m), and the Bamoun Plateau (1000 to 1300 m). Ages of the volcanic products along these edifices range from the Eocene to the Pliocene periods (Moundi et al., 2007;Moundi et al., 2008;Ngongang et al., 2015). The Noun and Ndop plains are flatted terrain, with average altitudes of 1000 m while the Mbo plain's altitudes range from 700 to 800 m, and there are several volcanic lakes. The geological formations made up of products of volcanic eruptions are lying on a basement rocks of plutono-metamorphic nature, and of Precambrian to Panafrican age (Njilah, 1991;Njonfang et al., 1998;Djouka et al., 2008). These basement rocks are usually associated with basic rocks (amphibolite and monzodiorite) and are masked in some places by a thick volcanic cover.
Annual rainfall increases from 100 mm to more than 3300 mm in the southern part around the city of Bafoussam, and 20 mm to more than 2400 mm when evolving to the northern part around Njimom due to the high altitude (Local agro-meteorology offices, Annual reports). Twelve months average temperature is between 26°and 28°Celsius. The vegetation mixes highlands forest and sub-tropical savannah, depending on the rainfall and the sun exposition. The population is 1,720,047 inhabitants, with a density of 124 inhabitants per km 2 .
In addition to the rainfall spatial distribution (Fig. 1d & e), human activities and settlements such as agriculture or buildings, mainly occupy slopes or shallows, exposing populations to natural hazards. For instance, on the 4th and 5th September 2018, terrain cracking followed by blocks slides damaged dozens of houses in the city of Foumban (IGMR-Penaye et al., 2018;Fig. 2a), caused the delocalization of hundreds of inhabitants. More recently, during the night of 27 to 28 October 2019, a long and huge rainfall of about 36 h triggered a rotational to translational landslide in Bafoussam (IGMR-Kankeu & Ntchantcho, 2019; Fig. 2b), the deadliest in that area with 45 deaths, dozens of missing people and at least 100 houses destroyed. Moreover, since the 1950's, more than 136 deaths were recorded in the area (Tchindjang, 2012(Tchindjang, & 2013. This context justifies the present research, to support Cameroon's government in anticipating such events and planning efficient mitigation actions.

Methodology
Original steps of SLIP/DRIP methodology In 2019, National Aeronautics and Space Administration's (NASA) developers proposed the SLIP and DRIP methodology, to automate rainfall-induced landslide identification in Nepal, by using open-source imagery and without the use of proprietary classification software (Fayne et al., 2019). It is a two-sided approach, that uses Landsat-8 imagery satellite imagery data to approximate visible landscape changes, and precipitation data for the landslide event's timing. Python 3 programming language and Google Earth Engine cloud environment support the computations, and the data analysis is based on a spectral band analysis and combined with ancillary field data.
Sudden Landslide Identification Product (SLIP) algorithm takes advantage of spectral properties of vegetation, slope, land-cover type, and soil moisture in biweekly (16 days) time steps to identify the affected area by a landslide right after the event, based on fresh bareearth exposure, and predict areas potentially exposed to upcoming events. The spectral red band is computed for Landsat-8, band 4, and computed as a percentage between the 10 composed recent images before the landslide and the most recent post-event images. Areas with at least a 40% increase in red reflectance are considered bare-earth, affected by the studied most recent landslide or exposed to landslide according to this criterion. To pursue, the soil moisture is assessed by adapting Normalized Multi-Band Drought Index, NMDI, of Wang and Qu (2007), to Landsat 8 spectroscopy. Basically, the NMDI monitors the soil and vegetation moisture using the following expression (1): where R 860nm , R 1640nm and R 2130nm represent the apparent at-sensor reflectance absorbed in the NIR and two SWIR wavelengths of the MODIS sensor measurements. However, the integration of the Landsat-8 band 6 that is closest to R 1640nm , gave poor results, and only bands 5 and 7 are used. The last step integrates a Digital Elevation Model's (DEM) slopes generating and thresholding.
The slope values are extracted in degrees and their intervals are classified as follows: gentle (0°-20°), fairly steep (20°-35°), steep (35°-45°), very steep (45°-60°), and extremely steep (60°-90°). All the values ≥20°are considered to be landslide-triggering. Further, noticing that a predominant triggering mechanism for landslides is rainfall (Petley et al. 2005), the Detecting Real-Time Increased Precipitation (DRIP) model leverage of NASA's Global Precipitation Measurement (GPM) was (re)built. It provides precipitation data with a more precise temporal window of occurrence for each potential event (Fayne et al., 2019). The DRIP algorithm identifies the likely timing of rainfall's peak, i.e., day and/ or hour of the day, that has triggered the studied landslide event, and that matches the SLIP detected areas, every 16 days. Windows of 24, 48 and 72 h are tested to obtain continuous rainfall data and integrate into the model.

Data acquisition and preprocessing
This experiment was conducted in a desktop script environment of licensed software, Erdas Imagine 2020 version 16.6.0.1366, ArcGisPro version 2.5 and XLStats 2020.1.64570. Twelve Landsat 8 satellite images were downloaded from the United States Geological Survey website, and the landslide of the 28 October 2019 was fixed as origin. These were then distributed such as 11 before and one after the event (Appendix 1). Due to the important cloud cover in the rainy season, and its effect on interrupting the 16-days temporal resolution necessary for the processing, all the best images available were collected, several from the dry season, i.e., December to March, for at least two images per year. Basically, free Landsat images are all level-1 products, and delivered as digital numbers (DNs). The bands used here are 2 to 7, namely bands blue to SWIR2, with a spatial resolution of 30 m, and the panchromatic band (Band8) was used to rescale the spatial resolution to 15 m (Table 1).
Applying the Cosine Solar TAUZ (COST) radiometric calibration model of Chavez (1996) to the stacked image, blue-SWIR2, the digital numbers were converted from atsensor radiance to top-of-atmosphere (TOA) reflectance via solar correction, and rescaled from 64-bit to unsigned 8-bit. Therefore, atmospheric corrections and haze reduction have helped to remove other noises and then approximate values of surface reflectance. The last step concerned the topographic correction that had addressed altitude artifacts. For the purpose of rainy season's land cover estimation, a classification map, change detection image and area expand function were applied (Appendix 2). The other entry is the Shuttle Radar Topography Mission (SRTM) image of the area with a spatial resolution of 30 m. It was also downloaded from the USGS website and preprocessed by using the void fill method to create a Digital Elevation Model (DEM) and reduce the errors of commission in flat areas where landslides are unlikely, such as riverbeds, which may have similar red reflectance and moisture characteristics (Jiménez-Perálvarez et al. 2011;Fayne et al., 2019). Its integration into the model helps defining the slopes threshold for landslide triggering.
Another entry concerns the precipitation data. These data were combined from three main sources. The Tropical Rainfall Measuring Mission (TRMM) (Braun, 2011), Tropical Applications of Meteorology using Satellite data (TAMSAT) (Maidment et al., 2014) and some local meteorology services.

Adapted SLIP algorithm
The first step is defining the soil exposure, i.e., the percentage of non-vegetated land. Fayne et al. (2019) proposed it as a rate of change in the red band reflectance between the current image before the landslide and a composed image of the 10 red bands of the images before the landslide. The formula is expressed as follows (2): Where Red current is the most recent Red band during or just after the landslide and is the 10 recent red bands just before the landslide. Then, the images should be collected for consecutives 16-days. In the present study, regarding the gap of almost 10 months in the same year between two Landsat 8 usable images, the percentage formula described above was leading to infinite values. Then the red difference was modified to an Inverse Normalized Difference Vegetation Index, INDVI, to assess the non-vegetated land. The INDVI is proposed as the spectral difference between the red and the NIR wavelength, such as (3): After an average of the INDVI was computed for the 10 oldest images, referring to the landslide of the 28 October 2019 in Bafoussam. Then, the average INDVI was subtracted from December 2019 INDVI and the resulting image was normalized in percentage to obtain the fraction of non-vegetated land (4): Where INDVIn stands for the normalized INDVI image, min and max are the minimum and the maximum of the INDVI. Values starting at 40% were selected as indicators for soil exposure to landslides as proposed by Fayne et al. (2019). A binned image was then created, with 0 for vegetated areas and 1 for non-vegetated areas, i.e., bare soil. Then, the vegetated class was expanded with factor 2 to approximate the land surface coverage in the rainy season, according to the classification statistics (Appendix 2).
Further, the soil moisture was assessed by using two indices. On one hand, the Modified Normalized Multi-Band Drought Index (MNMDI) (Fayne et al., 2019) was computed between the near infrared (NIR) and the shortwave infrared (SWIR2) bands (5): To confirm and complete the soil moisture information, the hydrothermal index composite was computed on the other hand. This index is used to enhance soils, rocks and minerals, as well as vegetation cover at a regional scale, based on a multiple ratios approach computation between the visible and infrared wavelengths (Pour, 2014). The concerned band ratios are SWIR1/SWIR (6/7), Red/Blue (4/2), and NIR/Red (5/4), while the result is a three principal components analysis image (Erdas Imagine, 2008 (2019) as significant to trigger landslide. After performing the hillshade processing to better highlight summits and valleys, the slope image was extracted in degrees. A binned image was coded, 0 for slopes less than 20°, and 0.5 for slopes at or above 20°. The three conditioned layers binned values are in Table 2.
The three layers were integrated using a simple weighted linear combination (SWLC) to map the areas where the conditions were met for sudden landslides. There are eight different values corresponding to the LHZ codes (Table 3). The SLIP process is described in Fig. 3 and the layers are represented in Fig. 4.
To validate the whole SLIP process, verification sites were chosen for each occurrence site, such as 2 in koutaba, 2 in Foumban and 5 in Bafoussam. The matching with affected sites defines accuracy.
Only remains the triggering factor identified as a long and huge rainfall. DRIP algorithm helps to model it.

Adapted DRIP algorithm
The DRIP tool is adapted as the rainfall intensity and threshold corresponding to the SLIP landslide mapping. Monthly precipitations of the west-Cameroon were computed between 1948 and 2017 for Africa, completed for years 2018 and 2019 (Table 2). According to the Tropical Applications of Meteorology using Satellite and ground-based observations (TAMSAT) data, especially its TRMM Multi-satellite Precipitation Analysis (TMPA) datasets component and mapping models, daily rainfall for the Cameroons' westhighlands were between 6 mm and more than 10 mm between 1983 and 2012 (Maidment et al., 2014). The annual highest rainfall period is between the second decade of June and the first decade of October, with at least 15 mm to more than 25 mm per day (Maidment et al., 2014;Dezfuli et al., 2017). The rainfall data collected by the local agro-meteorology offices (October 2019) assess the rainfall of 28 to 29 October in Bafoussam up to 81 mm, in about 36 h, before the landslide. This represents 22% of the 384 mm recorded for that month (Appendix 4). Ten groups of rainfall records were defined between June and October, that is 50 observations (Appendix 4). The month's selection is explained by the fact that all the landslides in Cameroon occurred in that 5 months interval, corresponding to the full rainy season ( Table 4).
The rainfall increases from June to September with highest records in August, and decreases in October, before stopping in November. The trends are the same for the number of rainy days, although some local differences can be barely noticed between two zones. In addition, the rainfall zonation was done from the lower (zone 1) to the higher (zone10) records. Samples of zones 1 and 8 illustrate these two statements for 2019 (Fig. 5).
The general trend gives an average rainfall of 2615 mm for 1948-2018 and 2573 mm in 2019, representing respectively 79% and 78% of the 3300 mm maximum annual rainfall.  August represents 25% of the 5 months and October records about 14% (Fig. 6). Moreover, the 81 mm of rainfall preceding the October 2019 landslide event represent 3.1% of the 5 months and almost four times the total daily rainfall, as mapped by Maidment et al. (2014) and Dezfuli et al. (2017). From the Table 4 above (including years 2012 to 2017), approximately 53 landslides have occurred in Cameroon between 1954 and 2019, and were all directly related to huge rainfall. As Fig. 7 illustrates, the highest month frequency corresponds to September with 24 occurrences (45%), followed by August with 12 occurrences (22%), October with 11 occurrences (21%), while June and July record the lowest same frequency of 3 occurrences (6%).
Further, their statistical distribution and variability were assessed by using the Weibull distribution. This method calculates a cumulative distribution function (CDF) or a probability density function (PDF) using the following eq. (6): Where the parameters are γ the shape, μ the location and α the scale. Because the Weibull model studies strength and failure of a system in relation with time (Klein, 2009), this study assesses the stronger relationship to the failures between rainfall and days of rainfall for the period 1948-2018 average, and then for the year 2019, setting =0.5 . The rainfall and days of rainfall data were rescaled by using the following ratios: a) [monthly rainfall ÷ monthly total of rainy days] b) [monthly rainy days ÷ monthly total days] The probabilities of failure and success were defined for both sets of data such as (7, 8): With PF and PS representing the probability of strength to failure and the probability of success. Then, their probabilities of strength to failure are suitable to be used as the z value in the standard PDF computation, which is defined as a normal distribution. Purposely and based on the data, the averages (≈3 and ≈4) and the standard deviations (≈1) of the mean daily rainfall and days of rainfall were computed for the two periods, and the PDF curves were superimposed on each other to find intersections ( Fig. 8a & b).
From the intersection of the two PDF curves, a set of thresholds were defined in a conditioning algorithm, with six explanatory variables (X), i.e., before landslide  and one dependent variable (Y), i.e., after landslide (2019). The whole conditioning algorithm is written in eq. 9 (i-vii) as follows:  & Then the probability of rainfall-triggered landslide is calculated by performing the logistic function. This is a classification algorithm useful for predicting binary outcome (1/0, Yes/No, True/False) given a set of predictor variables. It allows computing a multivariate regression between a binary dependent variable and several independent variables (Atkinson and Massari, 1998). Multiple logistic regression assumes that observations are independent and the natural log of the odds ratio and the measurement variables have a linear relationship. The quantitative relationship between the occurrence and its dependency on several variables can be expressed in the form of a logit function (10): Where P Ev , is the probability of an event occurring. In this case, the event is the daily rainfall threshold to trigger the landslide, and then being equal to 1 in the interval (0, 1). Z Ev. is the linear relationship of the event's occurrence with independent variables, and it is expressed as (11): Where b 0 is the intercept of the model, the b i (i = 0, 1, 2, …, n) are the slope coefficients of the logistic regression model, and the x i (i = 0, 1, 2, …, n) are the explanatory variables. Figure 9 synthesizes the processing.
The validation of DRIP model was conducted through a double accuracy assessment. The first is the confusion matrix, extracted from the logistic regression's classification table, based on the true positive rate, TPR, defining the ratio of all positive cases correctly predicted, and the false positive rate, FPR, expressing the ratio of all negative cases that are incorrectly predicted to be positive, under a defined threshold value. They are formulated as (12 & 13): The second is the positive predictive value, PPV, and the negative predictive value, NPV, that are respectively Modelling rainfall-triggered landslides from the survival analysis perspective The method of Caine (1980) was used to model the relation of rainfall to landslides. This process suggests a general threshold that works for time paces between 10 min and 10 days, using the rainfall intensity (I, mm/hr) and duration (D, hr). Here, the quantity of monthly rainfall was expressed as a function of the duration (Q, mm/hr) in the following eq. (16): The values of α and β are defined by using the Cobb-Douglas regression model in the formulation (17): The second step consists in assessing the spatial variability of the phenomenon. A spatial autoregressive model (SAR) enables to decompose the spatial process for a known site s, based on a random variable Z (s) , as follows (18): Where a and b are similarly extracted as for α and β. Therefore (19): Where, μ(.) is the spatial characterization and ε(s) is a centered random variable or error, α and a result from computing exponential of the intercept, and LHZ (s) integrates the four binary codes corresponding to at least two conditions met (Table 3) plus the fifth spatially closest code.
For an unknown site, S 0 , a predictionẑ s 0 of Z s 0 is interpolated using the observations Z s 1 , Z s 1 ,… Z s n through the kriging process, expressed such as (20): Then follows the Cox proportional hazards model. Originally, this model is used on the medical field to assess the probability that an individual will experience an event (for example, death) within a small-time interval, given that the individual has survived up to the beginning of the interval (Cox, 1972). The methodology mostly looks at the probability that given hazards, as the opposite phenomenon of hazards, may occur for a given actual occurrence. Although these models were not originally oriented spatial, they have been progressively integrated in geospatial analysis. Recent applications concerned fire hazard probabilities (Cyr et al., 2007) and factors of a space colonization (Baudains et al., 2015).
Theoretically, the hazard function for this case study is expressed such as (21): It gives the hazard function at time t for any unknown site S 0 with covariate vector that are the known sites of occurrence Z s . There are five known sites that are, one (1) in Koutaba, three in Foumban (3) and one in Bafoussam (1). There are also nine verification or unknown sites, distributed such as two in Koutaba, two in Bafoussam and

Landslide Hazard zonation (LHZ) SLIP-based map
The SLIP outcome is a map identifying the landslide zones of occurrence in eight classes (Fig. 10). The lowest class value, 0, stands for none condition met and the highest-class value, 2.25, corresponds to the landslide full conditions. The study subset is widely exposed to landslides hazard at different degrees.  ; 2%). The none condition class, 0, as well as the single conditions classes coded 0.5 (Slope inclination), 0.75 (Soil moisture) and 1 (Bare soil), cover 1276.29 km 2 , representing 32.5% of the study subset (Fig. 11).

Triggering rainfall DRIP-based thresholds
The highest concentration of rainfall and days of rainfall reveals information related to landslides occurrence (Fig. 12).   The maps outputted from the rainfall data show how the three sites, Koutaba, Foumban and Bafoussam match the main rainfall hot spots. The largest and densified spot is in the southern area, the second one is in the west and the third one is in the northeast area. At first sight and according to a visual comparison with the SLIP algorithm output, it may be assumed that rainfall concentration, in terms of quantity and number of days, is related to landslides. The further analysis of the logistic regression helps to comfort and discuss this assumption. The prediction model of rainfall-triggered landslide thresholds is expressed as (22) As a reminder, X 1 represents the rainfall frequency and X 2 stands for the days of rainfall frequency between   1948 and 2018, while Y is their intersection for the year 2019 (Table 5).
There are six ranges of probability over nine (6/9) that are rainfall-triggered landslides event, Y. The landslide probability is obvious at 99.99% when 0.229 ≤ X 1 , X 2 < 0.394 or X 1 ,X 2 ≥ 0.394. However, the influence of the rainfall frequency on the landslide is higher than the days of rainfall frequency, such as for X 1 ≥ 0.229 and X 1 ≥ 0.394, the probability remains 99.99%, no matter if the rainfall frequency X 2 < 0.229. Reversely, for all rainfall frequency X 1 ≤ 0.229, the probability of landslide occurrence is very low, 4.88%, no matter if for the days of rainfall frequency, 0.229 < X 2 ≤ 0.394 or X 2 ≥ 0.394. Consequently, the frequency of rainfall alone is able to trigger a landslide event in the study area once the minimum threshold of 0.229 is reached. Therefore, the adapted DRIP approach shows suitability to distinguish landslide and no-landslide for one common frequency.

Lifetime of sites to rainfall-triggered landslide
Applying the Cobbs-Douglas formula, the following models were withdrawn from the analysis of variance, ANOVA, regressions on the five sites of landslides occurrence (Table 6).
These values were normalized by dividing each site S i value by the mean of all values for each period. Then after, the most recurrent interval between LHZ codes, 0.25, was used to increment or decrement the normalized values towards the closest map class. Further, the means of the originally normalized values (Table 7) were computed for each S i , giving 10 values for the study period 1948-2019. These means vary between 0.009 (Koutaba) and 3.7 (Njiloum2-Foumban) and served as coefficients for the kriging process to predict the verification sites S 0 (Table 7). In   addition, the nine verification sites S 0 were coupled in four dependent variables for the linear regression with the known sites S i , and their LHZ codes were ordered in 10 value enhancing the variability, such as: i)Foumban: verification sites 1 and 2= S 0 -Foumban; ii) Koutaba: verification sites 1 and 2= S 0 -Koutaba; iii) Bafoussam: verification 1 and 2 = S 0 -Bafoussam1, verification 3, 4 and 5= S 0 -Bafoussam2 (Table 8). Table 9 gives the four detailed kriging models for each S 0 as well as the results. With respectively 1.99 and 1.642, S 0 -Bafoussam1 and S 0 -Bafoussam2 are the highly exposed sites to a potential rainfalltriggered landslide. S 0 -Foumban records the third value, 1.06, while S 0 -Koutaba, shows the negative value − 8571. Further, based on the landslide occurrence dates, the following parameters were introduced in the Cox model computation (Table 10).
The start date corresponds to the first landslide archived in Cameroon. The end date is the landslide occurrence events on the site. The "point date" is based on the analyst observations of the phenomenon, and in this case, it was defined according to the peak of rainfall starting in August, as well as on the highest percentage of landslide's occurrence that are 22% in August, 45% in September and 21% in October ( Table 2). The elapsed time was estimated in days rounded to the upper bound unit. The "censored status" is 1 for the failure to survive, that is the occurrence after the "point date" and 0 for success, i.e., the success or no-failure to survive before the "point date".
From the results presented in Fig. 13a, the survival probability of a site to rainfall-triggered landslide under 23,300 days of age was ≈0.68 (68%), and the site S4-Koutaba is the only concerned in this category. Between approximately 23,370 days and 23,700 days of ages, this probability was ≈0.38 (38%). The three sites S1-Nji-loum1, S2-Njiloum2 and S3-Njitout in Foumban, belong to this interval. Above 23,700 days of age, the survival probability keeps decreasing between ≈0.38 and ≈0.1 (10%). The only site that matches this category is S5-Gouache in Bafoussam.
Inversely, and based on Fig. 13b, the rainfalltriggered landslides hazard increases with time. Technically then, the exposure of S5-Bafoussam to that hazard is higher than the other sites, while the lowest exposure is at S4-Koutaba. The hazard ratio is 1.474 ≈ 1.5. (Table 11), corresponding to the time-toevent, meaning that at any time, one-and-half as many sites of occurrence (i.e., 1.5 * 5 = 7.5) are exposed to landslides.
The results above commented and their beta coefficients, β = 0.388, were used to elaborate Cox   proportional hazards model for the four unknown sites, S 0 . Their hazard ratio instantaneous risk was computed using the models detailed in Table 12 below. The highest ratio is observed at S 0 -Bafoussam1 with 1.47 that is almost as same as the hazard ratio. The other sites in the decreasing order are S 0 -Bafoussam2, S 0 -Foumban and S 0 -Koutaba, with respectively 1.28, 1.02 and 0.49.

Accuracies assessment and caveats
The two algorithms of SLIP and DRIP such as originally developed (Fayne et al., 2019) and adapted in this research connect landslides occurrence and huge rainfall. The interest has been to adjust the process to gaps in data, precisely the discontinuity among satellite images due to cloud cover and the unavailability of rainfall data at an hourly pace.
According to the SLIP processing, the highly exposed class to landslide occurrence meets the three conditions of bare soil, soil moisture and slope inclination thresholds, i.e., the LHZ code 2.25. Based on their location, five landslide sites match or are closely surrounded by this class, confirming the efficiency of the mapping method. In addition, five out of nine (5/9) verification sites support this statement, and give with the five previous a total match up to 71% (10/14 sites). Nonetheless, these accuracies are affected by the satellite images season, because the areas calculated may be x time more or less referring to the appending processing applied to approximate vegetation versus no-vegetation area between the dry and the rainy season (Appendix 2). In addition, the built-up extent and material can introduce biases especially when computing soil spectral indices (Ngandam et al., 2019). For instance, cities as Bafoussam and Foumban are characterized by their mi-rural/mid-urban patterns that include many houses in raw material such as earthen bricks and straw roofs, or unpaved dusty/muddy roads and tracks. Therefore, their reflectance may create mixedpixel in the INDVI result, because they usually reflect enough in the red and SWIR wavelengths of Landsat 8 images just as the landslide-affected areas (Ngandam et al., 2019). Further, the DRIP modelling accuracy assessment holds in two approaches. Taking the confusion matrix extracted from the logistic modelling classification table, and with 50 observations of rainfall in the 10 zones, the TPR or sensitivity is up to 100% (39/39), while the FPR or specificity is 82% (Table 13). Both give an overall accuracy of 96%, corresponding to the rate of rainfall and days of rainfall frequencies equal or beyond thresholds triggering the landslides.
In addition, from the data of Table 13 above, the PPV is 95% (39/41) while reversely, the NPV is 5% (2/41). These are high accuracies of the efficiency of a post-landslide analysis in relation with a daily rainfall, to potentially know a date of occurrence. However, the unavailability of hourly rainfall to proceed to a timing limits to properly correlate the two data frequency.
On the last step dedicated to the survival analysis of sites, the goal was to predict and correlate daily rainfall, days of rainfall and the magnitude of the landslide in terms of speed of occurrence, to complete the status mapping and timely retrospective of the original algorithms. Based on the rainfall intensity deducted on a daily pace, the processing was able to perform a lifetime analysis, departing from the first event archived in Cameroon on 1/1/1954, for the known and unknown sites. The survival probability of affected sites decreases with time, while the hazard of rainfall-triggered landslide increases. The unknown sites Cox's proportional hazards model can then be applied to the other sites of Cameroon, where rainfall data are available and landslides historic are archived. However, the huge rainfall of 81 mm in approximately 36 h preceding the landslide of Bafoussam-Gouache was not especially integrated in the processing, what raises the interrogation on the accurate timing as well as the rainfall intensity threshold to be used in the model.

Conclusion
The goal of this paper was to propose a derivative SLIP/DRIP procedure, adapted to countries and world's areas facing issues such as, satellite images important and a permanent cloud cover, poor spatial  and temporal records of rainfall data for hourly or daily paces, or the lack of landslide inventories, all that impede the implementation of the robust SLIP/ DRIP algorithms. Clearly thus, the proposed methodology did not pretend to substitute the standardized SLIP and DRIP algorithms. This approach is just another entry to these two algorithms, as an alternative to data voids, with the complementary step of survival analysis. To sum up, the ante and post landslide status of the affected areas approximately mapped, the rainfall threshold that triggers the landslide were modeled/estimated, and warning thresholds for a potential upcoming landslide were predicted for unknown sites based on survival probability and/or hazard exposure of known ones. Nevertheless, and logically, the accuracies of the outcomes suffered drawbacks and caveats, related to the accuracy of the areas affected by landslide, an efficiency of monitoring exposed areas in a window of 16 days with free Landsat 8 images, and the exact timing of the landslide occurrence in relation with the triggering rainfall.

Appendix 2
Image classification and NDVI computation for LULC extent assessment Because all the cloud-free satellite images of the study area are only available for the dry season, there is a need of matching the land use land cover (LULC) areas with the rainy season when landslides always take place. Then, an image of the rainy season was used from 23 June 2019, and a cloud-free subset was extracted on the natural land covers area, i.e., vegetation and soil, and a supervised classification was performed for each image. Their average overall accuracy is 92% and the average kappa coefficient is 0.89. An image difference was then performed with the newest image (dry season) classification of the same year. The transformational technique that produces a change image from which a change/no Appendix 1 change threshold must be established was used. It is expressed as a symmetric relative difference in the following equation (Erdas, 2008): Where Va is the new vegetation area and T is the time season image. It was noticed that in June, the 337,845 ha of the classification subset were occupied by the vegetation up to 61% (205,295 ha) versus 39% (132,550 ha) for soils, while in December, these percentages switch to 46% (155,260 ha) for vegetation and 54% (182,587 ha) for soils. To confirm the objects extraction and the trends above, the Normalized Difference Vegetation Index (NDVI) was computed for the two images (Fig. 14).
Statistics give 69% (233,113 ha) for vegetation and 31% (104,732 ha) for soil in June, versus 47% (158,787 ha) for vegetation and 53% (179,058 ha) for soil in December (Fig. 15). The average percentages are 65% for vegetation and 35% for soils in June, versus 46.5% for vegetation and 53.5% for soils in December. The ratios of the rainy season over the dry season areas were computed, showing that the rainy season vegetation area is about 1.4 times bigger than in dry season. Assuming that the average percentages could have the same influence on the classification process, the ratio of vegetation extent (65%) over the classification accuracy (92%) was calculated. The result obtained, i.e., 0.598 ≈ 0.6, was summed with the previous value, 1.4, as the best vegetation extent approximation for the rainy season, i.e., 2 times the vegetation area of the dry season's area. The ArcGISPro software expand function tool is useful for this purpose. In its principle, the class value targeted is multiplied by an x factor (2 here) to approximate the area as needed.
The algorithm is written as: Expand in raster; number cells; zone values ð Þ With in_raster representing the reclassified raster image, number_cells being the x factor and zone_values standing for the class to be expanded.