*2.1. Study Area*

The Wusuli River Basin is situated in the cross-boundary zone of Northeast China and the Far East region of Russia. Based on the Shuttle Radar Topography Mission (SRTM) 90 m Digital Elevation Model (DEM) (http://srtm.csi.cgiar.org/index.asp), we defined the entire boundary of the Wusuli River Basin with spatial analyst tools in ArcGIS 10 [34] (Figure 1).

**Figure 1.** Location of the Wusuli River Basin.

The latitude of this area ranges from 43◦25N to 48◦56N and the longitude ranges from 129◦50E to 138◦05E. The watershed area is about 195.06 × 10<sup>3</sup> km<sup>2</sup> in total, of which the Russian and the Chinese sections account for 69.06% and 30.94%, respectively. According to the administrative boundary, the basin runs from west to east across the Heilongjiang province of China and the Primorsky Krai province of Russia. In terms of the natural landscape, the Chinese section of the Wusuli River Basin is located on the East Sanjiang Plain, and the Russian section is part of the West Sikhot Mountains. The average elevation of this basin is 354 m, and the climate is characterized by a cold and dry winter and a warm and rainy summer.

The Wusuli River forms the border between China and Russia in the Wusuli River Basin, and the shared boundary stretches about 901.34 km. In 2015, the Chinese population in the basin was approximately 4,213,763 (http://data.stats.gov.cn/english/), and the Russian population was 878,007 (http://www.gks.ru/). Agriculture and coal mining are the main industries in the Chinese section, and more than 50% of the region is a plain (i.e., Sanjiang Plain), which is one of the most vital grain production bases in China. More than 70% of land is covered by forests in the Russian section, and timber and mining are the major industries—agricultural land covers less than 6% [35].

In the basin, wetlands serve as a stopover and nesting area for substantial migratory and waterfowl bird populations, such as *Grus japonensis*, *Ciconia ciconia*, *Larus ridibundus*, *Aix galericulata*, and *Tetrao tetrix* [36]. In addition, they play a vital role in stabilizing regional water supplies, ameliorating floods and drought and purifying polluted water. Furthermore, fish harvesting from wetlands is a significant economic resource for regional communities. Therefore, wetlands in the Wusuli River Basin are of grea<sup>t</sup> value for ecological balance, sustainable development and human well-being.

#### *2.2. Data Preparation and Fieldwork*

Cloud-free Landsat Thematic Mapper (TM) and Operational Land Imager (OLI) images (30 m spatial resolution) were chosen as the basic data sources to analyze the temporal and spatial wetland dynamics for the Wusuli River Basin in 1990, 2000 and 2015. These images were obtained in the growing season from May to September to minimize the effect of seasonal variations on the accuracy of land cover classification. Each image was acquired for the same month or the same vegetation growing period between 1990 and 2015 (Table 1). All images were geo-rectified with the registration error being less than half a pixel and atmospherically corrected using the Fast Line-of-sight Atmospheric Analysis of Spectral Hypercubes (FLAASH) [37].


**Table 1.** Specific description of each image used for land cover classification.

Data for the annual average temperature and annual precipitation from 1990 to 2015 were collected at 92 meteorological stations in and around the study area (Figure 2). This allowed us to analyze climate change factors driving wetland change. These meteorological data were interpolated to obtain a spatially continuous surface. The choice of spatial interpolation methods is referenced in Lu et al. [38].

**Figure 2.** The location of meteorological stations, field survey points and ground reference data sites.

From 2012 to 2015, 315 ground survey points were collected in the Chinese section of the basin watershed. Owing to limited accessibility to the Russian section of the basin, visual inspection of high-resolution images from Google Earth, online photos and literature searches were carried out to collect land cover information during the period 2010–2015. This generated an additional 256 reference points. All field survey points and the reference data sites were used to evaluate the accuracy of the land cover classification results in 2015 (Figure 2). Owing to the lack of field survey data in 1990 and 2000, 600 independent points for each year were created by spatial analysis of create random points in ArcGIS 10 [34]. These random points were classified into different land cover types (described in Section 2.3) by consulting with experienced interpreters and experts, and were then used as validation points.

#### *2.3. Land Cover Classification System*

Considering systematically the wetland classification of Ramsar Wetland Convention, the purpose of our study, and the specific conditions of the land cover type in the study area, a landscape classification system was established for this study, including nine land cover types (i.e., swamp/marsh, natural open water, human-made wetland, woodland, grassland, paddy field, dry farmland, built-up land

and barren land). These were incorporated into seven categories (i.e., natural wetland, human-made wetland, woodland, grassland, cropland, built-up land and barren land). The detailed description of the classification system is given in Table 2.


**Table 2.** Description of the landscape classification system used by this study.

#### *2.4. Rule-Based Object-Oriented Classification Method*

Compared with the frequent generation of 'salt-and-pepper' effects based on pixel-based classification methods [39], the object-oriented classification method can not only effectively avoid the 'salt-and-pepper' effects, but also reduce the 'within-class' spectral variation through segmenting an image into groups of contiguous and homogeneous pixels (image objects) as the mapping unit [40]. Moreover, besides of the spectral properties of the objects, their shape, texture and geometric features are also taken into account in the classification process of the object-oriented classification [8,41]. As a result,moreeffectiveandaccurateperformancesareobtainedthanwithpixel-basedapproaches[42,43].

To develop land cover maps for the study area in 1990, 2000 and 2015, a rule-based object-oriented classification method was applied to perform image segmentation and classify image objects into specific land cover types in the study. The layers which were selected to segmen<sup>t</sup> are Band 1 (0.45–0.52 μm), 2 (0.52–0.60 μm), 3 (0.63–0.69 μm), 4 (0.76–0.90 μm), 5 (1.55–1.75 μm), and 7 (2.08–2.35 μm) for Landsat TM, as well as Band 2 (0.450–0.515 μm), Band 3 (0.525–0.600 μm), Band 4 (0.630–0.680 μm), Band 5 (0.845–0.885 μm), Band 6 (1.560–1.660 μm), and Band 7 (2.100–2.300 μm) for Landsat OLI. The eCognition Developer 8.7.1 [44] was used to classify the images.

First, an optimal segmentation scale model referenced by Lu et al. [45] was used, in which a selected image scene was processed and grouped into homogeneous pixels (image objects) with an optimal segmentation scale. Each object resulting from this segmentation had minimal spectral variability [40,46] and the boundaries of these objects approximately followed the outline of individual land cover types. After segmentation, the segmented objects were categorized using a set of classification rules.

Considering the importance of vegetation growth and water content in wetland classification [47], the normalized difference vegetation index (NDVI) [Equation (1)] and land surface water index (LSWI) [Equation (2)] were used as rule layers to characterize vegetation and background soil, respectively.

$$NDVI = \frac{\rho\_{air} - \rho\_{rad}}{\rho\_{air} + \rho\_{rad}} \tag{1}$$

$$LSWI = \frac{\rho\_{\text{nir}} - \rho\_{\text{suir}}}{\rho\_{\text{nir}} + \rho\_{\text{suir}}} \tag{2}$$

where ρ*red*, ρ*nir* and ρ*swir* are the reflectance values of Landsat TM Bands 3, 4 and 5 and Landsat OLI Bands 4, 5 and 6, respectively.

Several previous studies have reported that the hue of di fferent band combinations can be a crucial factor for identifying di fferent land cover types [48]. In this study, the hue was derived from a combination of Landsat TM Bands 5, 4, and 3 or Landsat OLI Bands 6, 5, and 4. The value range of the hue is between 0 and 1.

Di fferent land cover types present distinct image textures, which is another variable necessary for land cover classification [49]. Based on the Haralick algorithm and gray level co-occurrence matrix (GLCM), the texture homogeneity ranging from 0 to 1 of each object was calculated [44]. The higher the value is, the higher the homogeneity.

The shape of an image object is also important to detect di fferent land cover types. The shape index (SI) [Equation (3)] of an image object describes the smoothness of an image object border. The smoother the border of an image object is, the lower its shape index.

$$SI = \frac{b\_v}{\sqrt[4]{P\_v}}\tag{3}$$

where *bv* is the border length of each image object, and *Pv* is the area of each image object.

After a series of pre-experiments, a classification rule set was developed (Figure 3). When the execution of classification rules was completed, the results were visually examined and modified for better precision. The overall accuracy, user accuracy and producer accuracy were used to assess the accuracy of the classification results.

**Figure 3.** Rules for the land cove classification of the Wusuli River Basin (<sup>α</sup>, β, γ, δ, ε, and ζ, represent the selected classification parameters and each of them could vary for different images (Landsat TM/OLI: Landsat Thematic Mapper/Operational Land Imager; NDVI: Normalized di fference vegetation index; LSWI: Land surface water index; SI: Shape index; SWIR1: Short wave infrared band1, corresponding to Landsat TM Band 5 and Landsat OLI Band 6, respectively; NIR: Near infrared band, corresponding to Landsat TM Band 4 and Landsat OLI Band 5, respectively; Red: Red band, corresponding to Landsat TM Band 3 and Landsat OLI Band 4, respectively; Unwater area: Land areas covered without water, including built-up area, barren land, woodland, grassland, swamp/marsh, dry farmland and paddy field).

#### *2.5. Analysis of Land Cover Change*

Two indices, annual land change area (ALCA) and annual land change rate (ALCR), were calculated to assess the dynamic degree of land cover types objectively. These are defined as follows:

$$ALCA = (\mathbb{U}l\_b - \mathbb{U}l\_a) \times \frac{1}{T} \tag{4}$$

$$ALCR = \frac{Ul\_b - lI\_d}{lI\_d} \times \frac{1}{T} \times 100\% \tag{5}$$

where *Ua* and *Ub* represent the area of each land cover type at the beginning and the end of the study period, respectively, and *T* is the number of years. In the study, the time interval was divided into two stages: 1990–2000 and 2000–2015.

To analyze the spatial change characteristics of natural wetlands (i.e., swamp/marsh and natural open water) more explicitly, intersect overlay analysis in ArcGIS 10 [34] was used to create a conversion matrix between natural wetlands and other land cover types for the time periods 1990–2000 and 2000–2015. In addition, a Sankey diagram [50] was used to illustrate the conversion of all land cover types, as this can help visualize the temporal dynamics of all land cover types.

#### *2.6. Calculation of Landscape Metrics*

Landscape metrics can reflect the characteristics of changing landscape patterns, and allowed us to assess quantitatively the landscape change process. In the study, five landscape metrics were used to assess the change pattern of the natural wetland landscape, including the number of patches (NP), mean patch size (MPS), largest patch index (LPI), area-weighed mean shape index (AWMSI), and the interspersion and juxtaposition index (IJI).

NP is defined as the count of patches and is a simple measure of fragmentation of one landscape category. Although the NP of one landscape category may be important for ecological processes and landscape pattern, it cannot directly reflect information concerning the distribution, area and density of patches. MPS is defined as the average patch size, and LPI quantifies the percentage of the largest patch accounting for the total area of all patches belonging to a given landscape category. AWMSI is used to assess shape characteristics by calculating the sum of the area-weighted ratio between the perimeter and area of each patch. IJI represents interspersion and juxtaposition and can quantify the connectivity and distribution pattern between di fferent patch types.

The detailed ecological significance and equations for selected landscape metrics are illustrated in Table 3. The calculation of landscape metrics was performed in Fragstats 4.2 [51].

#### *2.7. Climate Change Analysis Based on Mann–Kendall Test*

To measure the possible influence of climate change on the existence of wetlands, the changing climate trends were analyzed to determine whether climate change a ffected wetland dynamics. The statistical significance of the trends in annual average temperature and annual precipitation was measured using the Mann–Kendall test [52]. A trend is statistically significant if it is significant at the 5% level.



type; *aij* is the area of the *ij*th patch; *Pij* is the perimeter of the *ij*th patch; *eij* is the total length of edges between the *i*th and *j*th land cover type.
