## 1 Introduction  
In recent years, the tension between climate extremes and urban development has become increasingly prominent [1]. Flooding caused by extreme precipitation not only threatens the safety of people's lives and properties but the large amount of precipitation may also lead to soil saturation and foundation settlement, which poses a safety hazard to buildings [2]. Accurate and timely updating of the hazardous area or susceptibility maps is an important guide to master disaster data and enhance the ability of emergency response to prevent and mitigate disasters. Remote sensing technology has enabled the implementation of multi-scale and multi-frequency monitoring methods: from small-scale detection at high update frequencies to large-scale monitoring at low update frequencies and providing high update frequencies when necessary. Satellite-based synthetic aperture radar interferometry (InSAR) is an effective satellite-based remote sensing method for monitoring surface deformations over large areas with millimeter- to centimeter-scale vertical accuracy and bi-weekly or monthly temporal resolution [3]. However, the application of D-InSAR in monitoring surface deformation is usually affected by incoherence and atmospheric effects. In this context, multi-temporal InSAR processing methods have been proposed, such as the permanent scatterer interferometry (PS-InSAR) method [4] and small baseline set (SBAS-InSAR)methods [5]. For urban surface deformation, PS-InSAR is one of the most widely used and mature techniques in the field [6]. For geologically stable cities, surface deformation is a slow process of change, and potential areas of deformation hazards are difficult to identify in a short period of time using existing hazard assessment indicators. For example, the deformation of a slowly subsiding site tends to show a linear downward trend with a certain periodicity in a short inter-annual time period, but when the surface deformation values of a site repeatedly show anomalies beyond the original evolutionary characteristics, the site can be considered as a surface deformation-sensitive area. For example, [PERSON] et al. [7] identified deformation anomalies by counting the residual terms decomposed from a dam deformation time series. [PERSON] et al. [8] used a deep learning method to examine the deviation of the expected value of surface deformation from the actual value and judged values larger than the standard deviation of all deviations from the statistics as anomalies. Although surface deformation anomalies may not vary significantly in numerical terms, increasing the number of anomalies can lead to structural fatigue and pose safety hazards for buildings. Therefore, for a city like Zhengzhou, which has been geologically structurally safe and persistently stabilized by human impacts in recent years, surface deformation anomalies can be used to identify hazardous or susceptible areas [9]. In statistics, an outlier is a data point that does not belong to a particular group; it is an anomalous observation that is very different from the other values and departs from a well-formed data set. There are various methods for outlier testing, and the mainstream ones are as follows: 3 sigma, Z-score, boxplot, Grubbs hypothesis testing, KNN, LoF, COFsos, DBSCAN, iForest, PCAA, utoEncoder and so on. Each of these methods has its own advantages and disadvantages, among which the boxplot method has been widely used in the study of anomaly testing by virtue of its computational simplicity, high efficiency and the reduced influence of anomalies (excellent detection of anomalies on sea surface temperature, spatio-temporal precipitation and general circulation model data) [10]. The robust skewed boxplot method can be improved to be more accurate in detecting anomalies in skewed distributions [11]. However, if statistical methods are to be used to identify outliers in time series, it is necessary to ensure that the time series is stationary.