LIANG Peng ,QIN Cheng-zhi , ,ZHU A-xing , , ,HOU Zhi-wei ,FAN Nai-qing ,WANG Yi-jie

1 State Key Laboratory of Resources and Environmental Information System,Institute of Geographic Sciences and Natural Resources Research,Chinese Academy of Sciences,Beijing 100101,P.R.China

2 University of Chinese Academy of Sciences,Beijing 100049,P.R.China

3 Jiangsu Center for Collaborative Innovation in Geographical Information Resource Development and Application,School of Geography,Nanjing Normal University,Nanjing 210097,P.R.China

4 Key Laboratory of Virtual Geographic Environment,Ministry of Education,Nanjing Normal University,Nanjing 210023,P.R.China

5 Department of Geography,University of Wisconsin-Madison,Madison,WI 53706,USA

Abstract Selecting a proper set of covariates is one of the most important factors that influence the accuracy of digital soil mapping(DSM).The statistical or machine learning methods for selecting DSM covariates are not available for those situations with limited samples.To solve the problem,this paper proposed a case-based method which could formalize the covariate selection knowledge contained in practical DSM applications.The proposed method trained Random Forest (RF) classifiers with DSM cases extracted from the practical DSM applications and then used the trained classifiers to determine whether each one potential covariate should be used in a new DSM application.In this study,we took topographic covariates as examples of covariates and extracted 191 DSM cases from 56 peer-reviewed journal articles to evaluate the performance of the proposed case-based method by Leave-One-Out cross validation.Compared with a novices’ commonly-used way of selecting DSM covariates,the proposed case-based method improved more than 30% accuracy according to three quantitative evaluation indices (i.e.,recall,precision,and F1-score).The proposed method could be also applied to selecting the proper set of covariates for other similar geographical modeling domains,such as landslide susceptibility mapping,and species distribution modeling.

Keywords:digital soil mapping,covariates,case-based reasoning,Random Forest

1.Introduction

Digital soil mapping (DSM) is a widely-used way of predicting the spatial distribution of soil based on the soil–environment relationship (McBratneyet al.2003;Zhuet al.2010,2018;Zhanget al.2017).Selecting a proper set of covariates that could effectively present the soil–environment relationship is a key step (often the first) in DSM modeling.Currently,the selection of covariates is completely dependent on the user’s assignments according to their knowledge on DSM covariate selection,while an improper selection of covariates will impact the quality of DSM (Zhu 2006;Behrenset al.2010;Luet al.2019;Huanget al.2020).Meanwhile,there exist a large number of potential covariates without clear guidance on which condition each potential covariate should be selected for a specific DSM application.It is hard for users (especially the non-expert users) to be familiar with every potential covariate and to assess their feasibility for the specific DSM application at hand.This issue is more prominent when the users of DSM are increasingly extended to include those with less DSM knowledge,such as the hydrologist,ecologist,and natural resources managers and planners (Rossiteret al.2015;Jianget al.2016;Dharumarajanet al.2019).Therefore,automatic determination of a proper covariate set could not only enable the non-experts to use DSM more easily,but also improve the work efficiency of DSM users.

Currently,statistical or machine learning methods are the main solution to assist users in the selection of covariates for DSM.These methods can be grouped into three main ways,i.e.,filter,wrapper,and the embedding (Guyon and Elisseeff 2003;Chandrashekar and Sahin 2014).In the filter way,users analyze the relationship between each potential covariate and the soil properties firstly,then select the covariates closely related to the soil properties for DSM.This way commonly uses filter methods such as Pearson’s correlation analysis (Lagacherieet al.2013;Wanget al.2018),and moment correlation analysis (de Carvalho Junioret al.2014).In the wrapper way,different subsets of potential covariates are created step-by-step (e.g.,progressively adding new covariates,or gradually eliminates the insignificant covariates) and tested with a specified wrapper method (e.g.,the generalized linear models),respectively,then the subset that maximizes the performance of the wrapper method will be selected.The methods used in this way include stepwise regression procedure (Zhuet al.2015)and recursive feature elimination (Shiet al.2018) and so on.The wrapper way selects the covariates based on the final performance of the wrapper methods,and therefore usually has a better performance than the filter way (Zhouet al.2016).However,the methods used in the wrapper way should be trained for each subset of the potential covariates,which costs lots of training time.Furthermore,the order of the subset submitted into the wrapper method directly influences the results of the wrapper way (Derksen and Keselman 1992;Larket al.2007).To overcome the shortcomings of the wrapper way,the embedding way is proposed to incorporate the covariate selection as one part within the process of training statistical or machine learning models,i.e.,using statistical or machine learning models to select the final covariates from the potential covariates directly.As the representative methods of the embedding way,decision trees (Greveet al.2012),cubist (Adhikariet al.2014;Maet al.2017),and Random Forest (Vaysse and Lagacherie 2015;Ramcharanet al.2018) have been used to select DSM covariates.Recently,the potential of integrating diverse statistical or machine learning methods has been also explored,e.g.,Luet al.(2019) integrated Pearson correlation analysis,generalized additive models,and Random Forest for selecting DSM covariates.

All the above mentioned statistical or machine learning methods are data-driven methods which need users to collect as many soil samples as possible.However,these methods usually fail when the available samples are too limited to effectively fit the statistical or machine learning methods (Zhuet al.2014).

There exist lots of practical DSM applications in which DSM experts built machine-understandable DSM models(including model selection,covariate selection,data preparation,and so on) based on their relevant knowledge and the application context (such as the soil property to be mapped,the characteristics of the study area,etc.).Taking advantage of the covariate selection knowledge within those practical applications could be one possible way to solve the problem of selecting DSM covariates when there exist limited soil samples.Thus,this paper proposes a case-based method to formalize the covariate selection knowledge contained in the practical applications by cases,which could automatically select proper covariates for the new DSM applications.

2.Materials and methods

2.1.Basic idea

The covariate selection knowledge contained in the DSM application is tacit and non-systematic,which means that such knowledge is difficult to be formalized into explicit rules or mathematical equations.

Cases serve as a suitable way to formalize prior,nonsystematic knowledge in the artificial intelligence domain(Kasteret al.2005).A case has two parts,i.e.,the problem part that describes application context information (Qinet al.2016),and a corresponding solution part.The common case-based reasoning method is to find previous cases which are similar to the current application,and then to take the solutions of those similar cases (often that of the most similar cases) as the solution for the new application.A case-based method of determining proper parameter values for digital terrain analysis has been proposed by Qinet al.(2016).Their primary evaluation results showed that the case-based method is a promising way to apply existing application cases to the task of accounting for application context and relieve some of the user’s burden during digital terrain analysis modeling.

Therefore,a case-based method could be developed to formalize the covariate selection knowledge in existing applications into the form of cases,then uses those cases to select covariates for new DSM applications.In such a way,non-systematic and tacit covariate selection knowledge could be used in an automated manner.According to the above basic idea,the proposed case-based method is designed in two main steps:case formalization,and casebased reasoning,described as follows.

2.2.Case formalization

The purpose of case formalization is to transform knowledge on covariate selection into a set of discrete cases.As mentioned above,each case has two parts,i.e.,the problem part and the solution part.The key to formalizing the problem portion of cases is collecting a set of factors that can effectively describe the application context information,and then designing a set of quantitative attributes to represent those factors which can be directly used in case-based reasoning.Formalizing case solutions,on the other hand,is to record the solution which corresponds to the particular case problem.

According to thescorpanframework proposed by McBratneyet al.(2003),the covariates used in DSM can be classified into seven types,i.e.,soil,climate,organisms,topography,parent material,age,and space.Of these,topographic covariates are the most widely used in DSM applications (McBratneyet al.2003).Topographic covariates can be derived from the digital elevation model(DEM) with GIS Software,and DEM is an increasingly easy-to-access data with a wide range of resolutions.When selecting topographic covariates,one needs to consider little beyond the characteristics of the study area,such as further considering the availability issue of the covariate dataset.With the increasing availability of different GIS tools,more and more topographic attributes have been successfully used as DSM covariates (McBratneyet al.2003;Bishop and Minasny 2016).However,because of the complexity of digital terrain modeling (Hengl and Reuter 2008) and the lack of guidance on the condition of using each individual topographic attribute as DSM covariate,non-expert DSM users who do not understand the functionality of each of the large number of available topographic attributes may select an arbitrary and sub-optimal set of topographic covariates for their DSM applications (Lecourset al.2017).Therefore,the current design of the proposed method is focused on the recommendation of topographic covariates,which illustrates the concept simply and without loss of generality.

Formalization of case problemsAccording to previous DSM studies,there exist some key factors that influence the selection of covariates (Grimmet al.2008;Greveet al.2012;Liuet al.2013;Moslehet al.2016;Zhouet al.2016;Dharumarajanet al.2017;Pahlavan-Rad and Akbarimoghaddam 2018).These factors are divided into two groups in this study,i.e.,mapping target description factors,and mapping area characteristic factors.Mapping target description factors include the soil property being mapped (e.g.,soil organic carbon,the percentage of sand,etc.),the soil layer being mapped (e.g.,A-horizon,the depth of 20–40 cm,etc.),and the mapping resolution (e.g.,the resolution of 10 m for DSM).Mapping area characteristic factors describe the characteristics of the study area,such as the size of mapping area,the complexity of the mapping area (such as terrain,climate,and land use),and so on.

A set of quantitative attributes are designed to describe these factors (Table 1).The soil property of interest is described based on a soil property classification (i.e.,physical or chemical properties) as per Ditzler (2017).The mapping soil layer is described by the top and bottom depths,which refers to the depths (in the unit of centimeter) from the top and bottom of the soil layer to the surface,respectively.The factors of resolution and area size for DSM can be described by their natural quantitative attributes,i.e.,the resolution and area size respectively.The current method only considers the complexity of terrain;for now,it does not include other factors which describe the study area,such as land use,geology,and climate.This is because the current design of the proposed method is focused on selecting topographic covariates for DSM,and the complexity of terrain is usually used for representing the topographic characteristics (e.g.,plain or mountain) which could often influence the selection of topographic covariates (Lecourset al.2017).In this study,a set of quantitative attributes including total relief,standard deviation of elevation,and mean slope (Qinet al.2016) are adopted to describe the complexity of terrain in the mapping area.

Table 1 Factors and attributes used to formalize the problem part of a case in this study

Formalization of case solutionsThe solution part of each case is a record of which covariates are used for the corresponding problem.Note that in real DSM applications,there may be different forms for the same covariate across multiple cases (e.g.,“slope” as an alias of “slope gradient”,or slope gradient with different units -degree or tangent value).These should be rectified so that each covariate has a single,consistent descriptor across all cases.

2.3.Case-based reasoning

In this study,the point of the case-based reasoning is to determine whether a covariate should or should not be used for a new DSM application,which is a binary classification task.The proposed method adopts the Random Forest(RF) method for this binary classification task.RF is an ensemble machine learning method which combines a number of randomized classification and regression trees(Breiman 2001).RF has many advantages,such as noise resistance,nonlinearity,having no requirements for variable distribution,and working well for imbalanced data.RF has been successfully applied to many fields,including DSM,as a mapping method (Grimmet al.2008;Henglet al.2015;Shiet al.2018).

Within the context of this study,RF selects case samples from the case base by bootstrap sampling to build classification trees for each potential covariate one-by-one.The final classification result (i.e.,whether the potential covariate should be used or not) is a vote by the trees in the forests for each covariate,weighted by their probability estimates.This means the classification result is the one with the highest mean probability estimate across the trees.

In the implementation of the proposed method,the RF was performed within the Python package scikit-learn.Scikit-learn is an efficient and open-source tool for data mining and data analysis (Pedregosaet al.2011).The use of RF in scikit-learn depends mainly on two parameters,i.e.,thentreeparameter,which sets the number of trees in the forests,and themtryparameter,which sets the number of features tried at each node.

2.4.Experimental data

The cases of existing DSM applications used in this experiment were extracted from peer-reviewed articles published in international scientific journals in DSM and DSM-related domains,includingGeoderma,European JournalofSoilScience,SoilScienceSocietyofAmerica Journal,Catena,GeodermaRegional,PlantandSoil,ScienceoftheTotalEnvironment,EcologicalIndicators,EnvironmentalMonitoringandAssessment,GIScience&RemoteSensing,andPLoSONE.Each article was manually selected according to three criteria.First,the article must list covariates used for a real DSM application.Second,the article must describe the extent of mapping area.Third,the article must completely state the mapping target,including the soil property to be mapped,mapping soil layer,and mapping resolution.

Each case from an article was formalized into case problem and solution parts,as per the design above.Each case problem was formalized by specifying the value for each attribute of the problem part.For some of these attributes,such as those describing mapping target,the values are easy to extract,because they are often explicitly declared in the article.However,the values of other attributes (mainly those of the mapping area characteristics,such as attributes for the terrain complexity) which were not explicitly declared in the article,need to be recalculated.To make the proposed method automatic and easy to combine with DSM tools,those implicit attribute values should be calculated by normal GIS processing methods or tools based on an existing dataset (such as free-access DEM,and so on).Considering that the SRTM DEM is almost the most frequently used free-access DEM in digital soil mapping,the attributes for describing the mapping area characteristics(i.e.,total relief,standard deviation of elevation,mean slope) were calculated based on the SRTM DEM with 90-m resolution for each case in this study,so as to keep the attribute values comparative among different cases.

The solution part of each case is the list of covariates(topographical covariates in this study) used in the DSM cases.Note that in most of the practical DSM applications,soil experts first selected a preset of potential covariates according to their domain knowledge,then used a statistical method to select those important covariates from the preset potential covariates,and finally adopted the important covariates for their DSM models (such as those based on RF and so on).In these situations,only these important covariates were recorded in the solution part of the corresponding case.Note that the same covariate may have different names in different articles,as mentioned above.For example,the catchment area may be called as contributing area,flow accumulation,CA,etc.,and the topographic wetness index may be called as TWI,topographic index,etc.There also may be different covariates which have highly consistent effects from the perspective of DSM (e.g.,topographic wetness index and SAGA Wetness index;aspect,aspect-north,and aspecteast;different algorithms for landscape classification).It is necessary to unify those alternative names for the same covariate,as well as those covariates with highly consistent effects,when the solution part of the case is recorded in the case base.For this paper,a total of 38 topographic covariates were recorded in the current case base.

In this experiment,a total of 191 cases were extracted from 56 articles (Fig.1).These digital soil mapping cases were conducted within 23 countries and regions,mainly including China (containing about 14% of cases),Brazil(14%),Australia (13%),France (13%),Iran (10%),Canada(6%),and India (5%).Table 2 summarizes the value range of each attribute used to describe the case problem parts in the case base used in this study.In general,the attribute value ranges in the case base are representative of the corresponding attributes of DSM applications.

2.5.The reference method for comparison

In order to evaluate the effectiveness of the proposed method,a method which simulates a commonly used way of selecting DSM covariates by novices (“novice method” for short) was used as the reference method.The assumption adopted in the novice method is that the more frequently a covariate is used in the case base,the more popular that covariate is in the DSM domain.Without considering the application scenarios of each covariate,the novice method just considers the using frequency of each covariate in the case base.

Fig.1 Spatial distribution of the cases used in this study (more information in Appendix A).

Table 2 Summary of the attributes describing the case problem part in the case base used in this study1)

Before selecting covariates for a new DSM application(i.e.,a validation case in this study),the covariates used in the case base will be sorted from more to less based on the using frequency of each covariate.Then,the novice method selects the most frequently used covariates in the case base,according to the number of covariates used in the original solution of the validation case.When several covariates have the same frequency in the case base,the first selected covariate should be the one that results in the highest performance (i.e.,the covariate actually used in the original solution of the validation case should be selected).

2.6.Cross-validation and the evaluation indices

The proposed case-based method was evaluated based on a cross-validation experiment using the Leave-One-Out(LOO) method.This means that there were 191 tests,such that each of the 191 cases in the case base was treated as a validation case (i.e.,a new application problem) in one test,while the remaining 190 cases were used as the training set.In each test,the corresponding training set is used to train a RF classifier for each of the 38 potential covariates under consideration,as described above.In each test,the trained RF classifiers were applied to the corresponding validation case,so to determine whether each of 38 covariates should be used for the new application problem.The resulting recommendation from each test was compared with the solution part of the corresponding validation case,and this comparison was used to evaluate the proposed case-based method.Due to the difficulty of collecting soil samples used in the peer-reviewed articles,the accuracy comparison between final DSM products using the predicted covariates and the original covariates respectively was not adopted to evaluate the performance of the proposed method in this study.

Three quantitative evaluation indices were used to evaluate the performance of the proposed method.

where the true positives (TP) is the count of covariates which were both recommended by the proposed method for a validation case and actually used according to the original solution for the validation case (in other words,the count of covariates correctly recommended by the proposed method);the false positives (FP) is the count of covariates which were recommended by the proposed method but not used in the original solution for the validation case;and the false negatives (FN) is the count of covariates which were actually used according to the original solution for the validation case but not recommended by the proposed method.Therecallindex represents the ratio of correctlyrecommended covariates from the proposed method to all actually-used covariates according to the original solution for the validation case.Theprecisionindex stands for the proportion of covariates correctly selected by the proposed method in all covariates recommended by the proposed method.TheF1-scoreis the harmonic average of the precision and recall indices.TheF1-scoreis within the value range from 0 (i.e.,the worst performance) to 1 (i.e.,the best performance with perfect precision and recall).

These three indices were also calculated for each individual covariate in order to evaluate the performance of the proposed method at the individual covariate level.To show the difference,they were namedrecallv,precisionv,andF1-scorev.For example,theprecisionvindex for an individual covariate is the ratio of validation cases in which the proposed method selected that covariate correctly to all validation cases where the proposed method recommended that covariate.

3.Results

Table 3 summarizes the statistics (i.e.,mean,median,minimum,maximum,and standard deviation) of each evaluation index for the results from the proposed method.Among the 191 tests of the cross-validation experiment,therecallvalues range from 0 to 1 with an average of 0.644 and a median of 0.667.This means the proposed method recommended 64.4% of the topographic covariates used in the original solutions,on average.For half of the validation cases,the proposed method performed with a recall value larger than 67%.The mean and median values of theprecisionindex from the proposed method are 0.704 and 1,respectively.For the validation cases coming from Santraet al.(2017),the proposed method correctly selected all topographic covariates used in the original solution for the validation cases and did not select any extraneous topographic covariates.The distributions of therecall,precision,andF1-scorevalues from the 191 tests of the proposed approach were also analyzed with a set of index value intervals,i.e.,[0,0.25],(0.25,0.5],(0.5,0.75],(0.75,0.9],and (0.9,1] (Fig.2).For each of the three evaluation indices,the count of the validation cases in the interval(0.9,1] is larger than that of any other interval under analysis.This shows that the proposed method performed effectively.

Table 3 compares the statistics of the evaluation indices for the results from the novice method with those from the proposed method.Note that therecall,precision,andF1-scoreindices are all the same for the novice method.This is because the novice method selects those covariates most frequently used in the case base for each validation case,according to the number of covariates used in the original solution of the validation case.This means that the FP value always equals to the FN value,thus therecall,precision,andF1-scoreindices from the novice method are all with thesame value.The mean and median values of the evaluation indices from the novice method are lower than those from the proposed method.In addition,for most validation cases,the novice method produced lower value evaluation indices than the proposed method (Fig.2).This shows that the proposed method performed better than the novice method.

Table 3 Statistics of the evaluation indices for the results from the proposed method and the novice method,respectively

Fig.2 Histograms of the evaluation indices for the results from the proposed method compared with those from the novice method.x-coordinate refers to the five intervals of the evaluation indices whose values range from 0 to 1,y-coordinate refers to the number of validation cases that located in each corresponding interval.

The performance of the proposed method for each covariate is further analyzed.Fig.3 shows the value ofrecallv,precisionvandF1-scorevindices respectively from the proposed method for each covariate.For all three evaluation indices,the proposed method reached a score larger than 0.6 for most of the covariates.In the case of theprecisionvindex,the proposed method reached scores larger than 0.833 for most of the covariates and performed with a meanprecisionvof 0.798 among all covariates (Table 4).This means that generally the proposed method performed well at the individual covariate level.

There are still several covariates for which the proposed method reached a lower value of therecallv,precisionvandF1-scorevindices (Fig.3).Most of these covariates(i.e.,slope height,exposition,and beers aspect) are rarely used in the case base.For example,beers aspect (which received a zero value at all three evaluation indices from the proposed method) was used by the original solution for just once case in the case base.When this case was used as the validation case in the cross-validation experiment,it was impossible to select the beers aspect as a covariate by the proposed method.

4.Discussion

4.1.Effect of imbalanced data

Will imbalanced data (i.e.,imbalanced use of covariates in the case base) affect the performance of RF classifiers as used in the proposed method? Khoshgoftaaret al.(2007)analyzed the performance of RF built from imbalanced data by experiments where the percentage of minority samples varied from 1.33% (highly imbalanced) to 30% (only slightly imbalanced).Their experimental results showed that RF is effective for imbalanced data.In this study we analyzed the distribution of therecallvindex from the proposed method among covariate groups divided by the times of each individual covariate used in original solutions in the case base,i.e.,[1,5),[5,10),[10,20),[20,40) and [40,121].As shown in Fig.4,the median of therecallvindex remains comparatively stable at higher values (i.e.,near 0.7) among the different covariate groups.This means that the proposed method with RF is effective for imbalanced data.Note that therecallvindex from the proposed method for covariates used 1–5 times in original solutions in the case base shows a high variability (0–1;see Fig.4).This is mainly because the numbers of times these covariates used in original solutions in the case base are so low that once one of these covariates was classified correctly or incorrectly,it would result in a very high or lowrecallvvalue,respectively.This situation is acceptable for these seldom-used covariates.

Fig.3 Value of the recallv,precisionv,and F1-scorev indices respectively from the proposed method for each covariate.Along the horizontal axis,the covariates are sorted by the times used in the case base in a descending order from left to right.DEM,digital elevation model;TWI,topographic wetness index;CA,catchment area;MrVBF,multi-resolution index of valley bottom flatness;ProfC,profile curvature;PlanC,plan curvature;DistVert,vertical overland flow distance;VD,valley depth;MS,mid-slope position;LC,landscape classification;DisAlt,altitude above channel network;MrRTF,Multi-resolution ridge-top flatness index;DistHorz,horizontal overland flow distance;DI,insolation;SPI,stream power index;TPI,topographic position index;DH,difference between the elevation from the 20-m DEM to 5-m DEM;FD,flow direction;SP,slope position;SL,slope length;TCI,transport capacity index;Distance,overland flow distance;TRI,terrain roughness index;LE,landform evolution;NH,normalized height;MBI,mass balance index;SH,slope height;TSI,Topo-slope index;CS,catchment slope;PI,protection index;CI,convergence index;BA,beers aspect.

Table 4 Statistics of the evaluation indices for each covariate from the proposed method

Fig.4 The boxplot of the recallv index from the proposed approach among covariate groups divided by the times of each individual covariates used in original solutions in the case base(i.e.,[1,5),[5,10),[10,20),[20,40),and [40,121]).Each box indicates the distribution (i.e.,median,minimum,maximum,and interquartile range) of recallv scores of the covariates within corresponding group of times used in original solutions.

4.2.Limitations of the proposed method

There still exist several limitations of the proposed method.Firstly,although the experimental results showed that the proposed method improved more than 30% accuracy than the reference method,the performance of the proposed method with a meanrecallof 0.644 in the experiment could be further improved.To improve the accuracy of the proposed method,a larger case base is necessary.Manual extraction of the cases in this study limits the size of the case base.Therefore,a semi-automatic or automatic case extraction method should be further developed to improve the efficiency of case base construction.Secondly,the detailed design of the case-based method proposed in this study focused on selecting topographic covariates for DSM.The proposed method should be extended to consider those factors and attributes that influence the selection of other types of DSM covariates (McBratneyet al.2003).

4.3.Potential application of the proposed method

While the application domain in this study was focused on DSM,the proposed method may also provide a new solution for other geographic feature mapping domains which often have a similar problem of selecting a proper set of covariates,such as landslide susceptibility mapping (Zhuet al.2014),species distribution mapping (Fourcadeet al.2018),and habitat suitability mapping (Harris and Baker 2011).This is valuable for developing easy-to-use geographical modeling tools and reducing the non-experts’ burden on both learning domain knowledge and complicated operations (Jianget al.2016;Qinet al.2016;Zhuet al.2020),especially when coupling with automatic data matching methods (Houet al.2019;Zhu and Yang2019) and automatic model determination methods (Jianget al.2019).

5.Conclusion

This paper proposes a case-based method for automatically selecting DSM covariates,which is based on the idea of using the tacit knowledge on covariates selection which was contained in existing practical DSM applications.Compared with the common statistical or machine learning methods for selecting covariates for DSM,the case-based method could effectively relieve the users’ burden of preparing lots of soil samples and a large number of pre-selected covariates.The experimental results show that the proposed method performed much better than the novice method.Thus,the proposed method is valuable for the users (especially the nonexpert users) on selecting covariates for DSM modeling.The proposed method could be also applied to selecting proper covariates for other similar geographical modeling domains.

Acknowledgements

This work was supported by grants from the National Natural Science Foundation of China (41431177 and 41871300),the Priority Academic Program Development of Jiangsu Higher Education Institutions (PAPD),China,the Innovation Project of State Key Laboratory of Resources and Environmental Information System (LREIS),China (O88RA20CYA),and the Outstanding Innovation Team in Colleges and Universities in Jiangsu Province,China.Supports to A-Xing Zhu through the Vilas Associate Award,the Hammel Faculty Fellow Award,and the Manasse Chair Professorship from the University of Wisconsin-Madison are greatly appreciated.