Effect of Land Use Change on Gully Erosion Density in the Black Soil Region of Northeast China From 1965 to 2015: A Case Study of the Kedong County

Kedong County is typical of the black soil region of northeast China in being highly susceptible to accelerated soil erosion by gullying. Using data sourced from Corona satellite imagery for 1965, SPOT5 for 2005 and GF-1 for 2015, the spatial distribution of gullies in the research area was mapped. Land use data for 1965, 2005, and 2015 were obtained from the topographic map of 1954, and from Landsat images for 2005 and 2015. Over the last 50 years, the extent of gully erosion in the study area has increased markedly, most notably on cultivated land, while gully density rose from 2,756.16 m2/km2 to 14,294.19 m2/km2. Cultivating land on slopes, especially on slopes greater than ∼4°, may rapidly aggravate gully erosion. The greatest increases in gully density occurred in situations when cultivated land and other/degraded land were transformed, which gully erosion density increased by 49,526.69 m2/km2. Other/degraded land is the most vulnerable land in the study area, with the highest gully erosion density. In these cases, gully density initially increases and, although the “Grain for Green” project has been implemented, gully erosion density has not always declined in the recent past.


INTRODUCTION
Over the next 50 years, the planet faces the daunting task of feeding 3.5 billion additional people (Borlaug, 2007). How can the associated increased demand for agricultural commodities be accommodated? Clearly, there will be very significant pressure on arable land globally and a need to increase its area (Gibbs et al., 2010;Amundson et al., 2015). Changing land use, however, presents us with a problem because it can dramatically alter the albedo and exacerbate climate change, accelerate rates of soil degradation and erosion, and cause loss of biodiversity. Soil erosion is one of the most prominent environmental impacts associated with land use change, and has resulted in extensive loss of productivity and ecosystem services at the global scale (Green et al., 2005;Montgomery, 2007;Ramankutty et al., 2008;Foley et al., 2011;Ercoli et al., 2020;Li et al., 2021). The associated decline in habitat quality and fragmentation caused by development has a wide range of deleterious impacts on soil conditions. It has long been recognized that erosion in excess of soil production eventually results in decreased agricultural potential (Liu and Diamond, 2005;Mukai, 2017;Macias-Fauria, 2018;Xiong et al., 2019).
China's agricultural enterprise is the largest in the world; the country produces 20% of the global grain harvest on 9% of the earth's arable and permanent cropland, which must profoundly influence global biogeochemical cycles (Frolking et al., 1999). However, soil erosion affects 19% of China's land area, one of the highest figures for any country (Liu and Diamond, 2005;Macias-Fauria, 2018). Among various types of soil erosion, gully erosion is most common and has both on-site and significant off-site ecological effects (Valentin et al., 2005;Cotler and Ortega-Larrocea, 2006;Zucca et al., 2006;Galang et al., 2007;Gutiérrez et al., 2009;Martínez-Casasnovas and García-Hernández, 2009;Dotterweich et al., 2012). Soil erosion has been reported to be especially prominent in the so-called "black soil" region of China (Valentin et al., 2005;Gu et al., 2018;Xu et al., 2018). The black soils of northeast China have been identified as being particularly susceptible to severe gully erosion (Xie et al., 2019). As the area under agriculture has expanded in this region, so has the frequency of gullies, with smaller rills developing into larger features. Indeed, Liu and Yan (2009) identified some 250,000 gullies here, accounting for a reduction in the area of cultivated land of some 4.83 × 10 5 ha, and equivalent to a loss of grain productivity amounting to 3.62 × 10 9 kg.
Gully erosion, as a sign of serious soil erosion and degradation, has been paid more and more attention (Amare et al., 2019). Soil degradation due to gully erosion is therefore an important work in ecological restoration efforts in the black soil region. Thus far, research on gully erosion in the black soil region has focused on three main aspects: viz., short-term, small scale monitoring of gully by differential GPS Wu et al., 2008;Hu et al., 2009); mapping the spatial distribution of gully erosion over time using remote sensing images, aerial photographs and GIS technology at relatively large spatial and temporal scales (Rysin et al., 2017) and dynamic monitoring of largescale (large area, long-term) and small-scale (small area, shortterm) gullies using multi-source, multi-scale remote sensing monitoring systems combined with ground measurements and high-resolution satellite, low-altitude UAV remote sensing (Pu et al., 2016;Feng and Li, 2018). The above studies have explored the influence of natural drivers, such as climate, topography, soil and vegetation, and human factors, such as choice of crop type, farming practices and shelter forest plantation (Deng et al., 2015;Yang et al., 2017;Wen et al., 2021). However, despite the importance of land use/cover change in aggravating soil erosion (Zhu et al., 2003;Rahmati et al., 2016;Mukai, 2017;Selkimaki and Gonzalez-Olabarria, 2017;Zabihi et al., 2018;Amare et al., 2019), there has been little research on the relationship between the development of gully erosion and land use change in this region (Li et al., 2016), although Wang et al. (2016) have done so recently at the national scale. Most of these studies consider the occurrence of gully erosion under different land use types, but few have attempted to examine the effects of land use conversion on gully erosion. For example, what type of land use conversion results in the greatest increase in gully erosion intensity? Do some types of land use change actually reduce gully erosion? The aim of this article, therefore, is to quantify the change of gully erosion since 1965, analyze gully erosion dynamics, explore the relationship of gully erosion to land use types, gully density dynamics and land use change, and establish the degree to which changes in gully erosion density are caused by land use change under different slope conditions at higher resolution. The study is applied to the Kedong County of Heilongjiang Province, which is considered as representative of the broader region and the results therefore reflect the importance of gully erosion in the black soil area of northeast China in general.

Study Area
Kedong County, Heilongjiang Province, is located between E126  (Figure 1). Kedong is in the cold temperate zone and experiences a semi-humid to semi-arid, continental monsoon climate. The climate in study area characterized by rainfall being concentrated in the summer months and by severe frosts in winter, which only permits cropping during the summer. The mean annual precipitation is 526.5 mm, with mean annual runoff of approximately 75 mm, derived largely from summer rains that fall mainly in July, August and September in the form of convectional storms associated with high rainfall intensity. The area is situated in the alluvial-diluvial undulating plain of the Xiaoxinganling mountain and Songnen Plain transitional zones. The main soils in the area are black soil, meadow soil and dark brow soil (which are Phaeozems, Haplic Phaeozems, and Haplic Luvisols in the WRB Classification System). The physical and chemical properties of the soils are uniform , and the erodibility factor (K) are 0.2512∼0.3108. Terrain is typical of the rolling hills of the black soil region (Fan et al., 2005). Elevation here ranges between 192 and 423 m above mean sea level, and most slopes, relatively long in nature, are less than 5 • . The landscape is vulnerable to accelerated soil loss, predominantly through gully erosion, due to a combination of factors including climate, soil and terrain. Agriculture is currently smallholder-based with widespread subsistence production. During the busy farming season, a mechanical tillage with a deep cultivation within 35 cm will be taken in a agricultural terraces of hillslopes. The main crops of the study area are soybean and corn, and the cultivated area is 58% and 41% of the total area, respectively. As an important commodity grain production area, the research area has long been subject to overuse and inadequate attention has been paid to the maintenance of soil fertility and structure (Wang et al., 2009;Duan et al., 2011;Fang and Sun, 2017;Gu et al., 2018). These circumstances have resulted in the deterioration of soil physical and chemical properties, loss of fertility, expansion of the area exposed to accelerated soil loss, and an associated increase in erosion intensity. Overall, there has been a reduction in grain production in the region which poses a threat to national and regional food security (Liu and Yan, 2009). Up to 75% of the land in this area has been converted to cultivated land (Yan et al., 2008) but inappropriate land use has led to accelerated soil erosion, a situation that has been recognized by the government and in response China initiated its 'Grain for Green' project in 2003.

Slope
The 1:50,000 topographic map in the study area was digitized to extract contour data and generate a digital elevation model (DEM) with a spatial resolution of 25 m by interpolation. The slope layer was generated from the DEM using the ArcGIS spatial analysis module. According to the slope classification of the black soil region defined in the Technical Standards for Comprehensive Control of Soil Erosion in the Black Soil Region (Ministry of Water Resources, China, 2009), the slopes were classified into six grades: <0.25, 0.25∼1.5, 1.5∼3, 3∼4, 4∼5, and >5, and a slope classification map of the study area was obtained (Figure 2). Spatial data including the digital elevation model (DEM), land use types and gully erosion were obtained from a range of sources (Table 1) for the time period 1965-2015.

Gully Erosion Data
The distribution of gullies was obtained by visual interpretation of remote sensing imagery, with spatial resolution as indicated in Table 1. Images were obtained from the CORONA  All images were taken in late June of each year, at the time when the cloud volume was less than 10%. Initially, the orthophoto maps from CORONA and SPOT5 were georeferenced against the 1:50,000 digital topographic map after control points were selected. The images from the GF-1 satellite were corrected using the referenced images from the SPOT5 imagery. In addition, to select control points based on the images from SPOT5, the control points were also ground-truthed through consultation with local informants. Following field investigation ( Figures 3A,B), establishment of interpretation keys, manual visual interpretation, field verification, re-interpretation and error amendment, distribution maps of gullies in the study areas for 1965, 2005, and 2015 were generated (Figures 4A-C). The spatial analysis of gullies in different periods shows that all gullies are active in the study area, and all of them expanded in past 50 years.

Land Use Data
The basic data sources for land use included the first generation 1:100,000 topographic maps compiled by the Surveying and Mapping Bureau of the Headquarters of the General Staff in 1954, and the land boundary was corrected by combining these with high-resolution images from CORONA for 1965. Land use data for 2005 and 2015 were obtained from Landsat/TM, and in each case for June when cloud cover is minimal. The data acquisition process comprised of, firstly, image correction using a mathematical simulation of geometric distortion on the original images using ground control point data, and establishing the degree of correspondence between the original image and the spatially corrected image. According to a slightly modified version of the LUCC classification system (Liu et al., 2014), six primary classes were identified, viz., cultivated land, forest land, wetland, water, construction land, and "other/degraded" (which is mainly abandoned, degraded or bare land, and may include   other vegetation types such as grassland). Based on the data from field investigation and expert knowledge, with reference to relevant geographic maps, the ground object spectra were analyzed to establish unified image interpretation keys. Visual interpretation was conducted by ArcGIS software to obtain the land use data in the study area (Figure 6).

Gully Erosion Dynamics
There were 767 gullies in the research area in 1965, with a total length of about 324 km; the area of eroded land was about 4.5 km 2 , and the maximum density was 68,316 m 2 /km 2 . In 2005, the number of gullies had increased to 2,322, with a total length of about 1,122 km and an area of eroded land of some 19.1 km 2 , with maximum density of 134,256 m 2 /km 2 . In 2015, the number of gullies had risen further to 2,746, with a total length of about 1,395 km; the area of eroded land was about 24.7 km 2 , and maximum density of 144,924 m 2 /km 2 . By 2015, the gullies had increased by a factor of 3.5 compared to 1965; equivalent ratios for length were ×4.3, for area ×5.5, and for maximum density approximately ×2. The increase gullies only contain single new gullies, does not include additional peaks from the main old gullies. Key: 1 cultivated land, 2 forestland, 3 wetland, 5 construction land, 6 other/degraded, 1→2 represents the conversion from cultivated land to forestland, and remaining likewise. According to Figures 5, 7, 1633.44 km 2 of the study area was subject to sporadic erosion in 1965, accounting for 78.71% of the total, while violent erosion covered about 38.56 km 2 , representing 1.86% of the research area. By 2005, the sporadic erosion area was reduced to 1030.11 km 2 , <50% of the study area; violent erosion now occupied 340.15 km 2 (16.39%), almost ten times that of 1965. By 2015, the sporadic erosion area was further reduced 899.98 km 2 (43.37%), and violent erosion reached 449.50 km 2 (21.66%). Overall, there was a substantial decrease in the area with no or low levels of gully erosion during the period studied, while, correspondingly, the area exhibiting extreme erosion increased markedly.

Land Use Change
We conducted statistical analysis of land use data in the study area (Figure 8) and established a land use transfer matrix. Table 2 illustrates land use transfers for the period  1965 to 2005, while Table 3 covers the 10 years period from 2005 to 2015. The main land use categories in 1965 included cultivated land, Wetland and forest land, accounting for more than 96% of the total area. Cultivated land was the dominant land use type and covered almost two thirds of the area, followed by wetland, which occupied around 20%, especially adjacent to the major rivers. Forest land accounted for about 10% of the total area, mainly distributed in the north and southeast. By 2015, cultivated land in the research area increased to 73.14% of the total area, while other/degraded land and construction land also increased significantly. Figure 8 indicates that the area of both forestland and wetland decreased markedly during the period under consideration. Table 2 shows the total area of land transferred between different land use types in the research area to be 500.84 km 2 between 1965 and 2005; approximately 25% of the total area was subject to land use change. Of the land use types, wetland underwent the largest change, having been reduced by more than 228 km 2 in total. The most common conversions were from cultivated land to construction land and degraded land.
The area of cultivated land, construction land and other/degraded in the region exhibited minor decreases between 2005 and 2015 (Table 3), while forest land and wetland increased, although overall land use structure did not change significantly. The total area transferred between different land use types during this period was 54.61 km 2 , accounting for less than 3% of the study area. Table 4 reveals that the gully erosion situation in the Kedong County has deteriorated over the period in question. In addition, it highlights that the area of cultivated land subject to erosion increased from 84.21% in 1965 to 87.81% in 2015. The area of eroded cultivated land increased by 17.90 km 2 , while gully density rose from 2,756.16 m 2 /km 2 to 14,294.19 m 2 /km 2 , representing an average annual increase of 230.76 m 2 /km 2 . Cultivated land was the land use type on which gully erosion increased most rapidly.

The Relationship Between Land Use Types and Gully Erosion
Gully density under the "other/degraded" land use category, increased from 3,390.52 m 2 /km 2 to 14,801.18 m 2 /km 2 over the 50 years period, with an average annual increase of 228.21 m 2 /km 2 ( Table 4). By 2015, the gully area in this land use category had increased by approximately by a factor of x18 relative to 1965, although there was a slight decline in between 2005 and 2015. The average annual rate of increase was just 56.33 m 2 /km 2 in forest land, which indicates that this land use type was generally effective at minimizing gully erosion. Wetland were less impacted by human activity, and gully area in this land use type declined markedly, especially during 1965 to 2005, this being mainly as a consequence of drainage and cultivation across substantial areas. Table 4 further shows that gully erosion in construction land continuously increased, at an average annual rate of 110.65 m 2 /km 2 during the past 50 years, suggesting that urban development aggravates gully erosion.

The Relationship Between Land Use Change and Gully Erosion Dynamics
By overlaying the land use change maps coupled with gully erosion distribution maps, the variation of gully erosion density related to land use type conversion can be obtained (Figure 9). Figure 9 reveals that it was the conversion of cultivated land and other/degraded that led to the most significant increase in gully erosion from 1965 to 2015. During this time 16.61 km 2 of other/degraded land was converted to cultivated land in the study area and gully erosion density increased by 49,526.69 m 2 /km 2 , with a mean annual change of 990.53 m 2 /km 2 . The area of cultivated land that was converted to other/degraded land (i.e., abandoned) was 30.18 km 2 , the density of gully erosion increased by 27,678.68 m 2 /km 2 , and the average annual change was 553.57 m 2 /km 2 . As shown in Figure 9, land use changes generally resulted in increased in gully erosion. For example, marked changes in gully density occurred following forest land and wetland being converted to cultivated land -15552.68 m 2 /km 2 and 12594.74 m 2 /km 2 , respectively -values which are 14x that of the original land types. The results show that forest land and wetland were easily eroded after being reclaimed into farmland, showing a trend of accelerated development of gully erosion. In addition, the increasing rate of gully erosion density failed to slow down even after cultivated land was converted to forest land, although the conversion of cultivated land to forest land would be expected to inhibit further development of gullies.

Gully Erosion and Land Use Change in Relation to Slope
By overlaying the land use change map, gully density map, and slope map, the variation of gully erosion density caused by land use change under different slopes during the period 1965 to 2015 can be obtained (Table 5 and Figure 10). Figure 10 shows that the conversion of forestland, grassland and wetland into cultivated land had increased the incidence of gullies, a pattern that is accentuated on slopes, especially those exceeding 4 • . Results of the statistical analysis indicate that, when cultivated land is changed to wetland, gully density increases only gradually on slopes of less than 3 • , but more rapidly when the slope angle exceeds 3 • ; This situation is similar when cultivated land is changed to other/degraded (abandoned). Generally, the conversion of other/degraded into cultivated land was the most significant cause of an increase in gully density on all slope. Gully erosion density increased when cultivated land changes is converted to forestland and wetland, especially on the steeper slopes (3 • ∼4 • ).

DISCUSSION
While land use/cover change is demonstrably an important factor in gully erosion (Arabameri et al., 2018,Amare et al., 2021, and affects the spatial distribution of erosion, research on how land use change impacts gully erosion intensity is less common. Following classification, validation and field verification of Landsat/TM images for 1965, 2005, and 2015, and this paper presents a detailed assessment of land use changes and its impact on soil erosion for a 50-year period in the black soil area in northeast China. The area is characterized by marked changes in both erosion intensity and gully length per km 2 during the last few decades and, indeed, demonstrates the highest rate of soil erosion in China (Liu et al., 2010;Gong et al., 2013;Wang et al., 2013). Compared with research results in other regions of the world (Amare et al., 2019), the emergence of gullies in this area is mainly located in cultivated land. Land use change associated with agricultural development in the black soil areas of northeast China has caused profound ecosystem change, and is the main driver of accelerated soil erosion in this region. Topography in the study area is characterized by gentle slope and long slope length. Due to improper farming on the low gentle slope farmland, gullies develop rapidly in the farmland (Wen et al., 2021). Through a comparison of gully length and area for 1965 and 2015; gully length in 2015 was increased by a factor of x4.3 over 1965, while the equivalent factor for gully area was x5.5. The rapid development of gully erosion has resulted in the fragmentation of cultivated land, which indirectly led to agricultural machinery cannot to cultivate within 3∼10 m on both sides of the gully to become abandoned area (Ou et al., 2018). The reduction in the area of cultivated land caused gully erosion is accompanied by the loss of soil organic carbon, a decline in soil quality and land productivity, and threatens food security (Garzon-Garcia et al., 2014;Basher et al., 2017). Once the soil fertility is reduced, independent farmers always abandon farmland and reclaim it in other places, which leads to initiation of gullies again. It can be seen that there is a positive feedback mechanism between land use change and gullies. In this study, we statistically analyzed the change of gully erosion intensity arising from changes in land use type from 1965 to 2015. The results showed that the gully erosion density change was greatest when cultivated land and other/degraded land were transformed, which also support that view.
Furthermore, the analysis reveals the vulnerability of other/degraded land, which often relates to farm abandonment and vegetation degradation; this has exhibited continuous increase over the past 50 years. When cultivated lands are abandoned and wetland degraded, changes in surface properties, especially the ensuant decrease in vegetation cover, leads to initiation of or renewal of gully erosion. This situation is particularly evident on slopes. There are few remaining natural grassland in the study area, and the scattered sparse grassland which is abandoned and degraded is especially prone to accelerated gully development. Some studies (see Gomez et al., 2003;Parkner et al., 2006;Wang et al., 2016) have shown that the catchment area required for gully formation in grassland is smaller than that for primary and secondary woodland, and that grassland is the land use type most vulnerable to erosion. Moreover, especially in the early stage of grassland conversion, soil permeability under the root layer is reduced, which raises the moisture content of the topsoil and renders it more easily prone to erosion.
The rapid increase in cultivated land and the reduction of forest land and wetland have all contributed to the intensification of gully erosion in the Kedong over the past 50 years, a situation which appears to be ongoing. The results of the first national water census showed that the developing gullies occur in up to 80% of the black soil area of northeast China, and there are approximately 295,633 erosion gullies of length >100 m in the black soil area of northeast China associated with an eroded area of some 3648 km 2 . This amounts to an annual loss of cultivated land of 7.39 km 2 annually; Although a series of environmental protection measures, such as "Grain for Green" , planting shelterbelts (Fang, 2017), mulching straw (Wen et al., 2015), terracing and contour tillage (Fang and Sun, 2017;Yan et al., 2021) appears to have arrested the development of gully erosion, these measures have thus far not been sufficient to control soil erosion comprehensively in the black soil region.

CONCLUSION
This article takes the Kedong County as an example to study the effects of gully erosion dynamics and land use change on the development of erosion gullies in a typical black soil area of northeast China during the past 50 years. The results show that: (i) In the past 50 years, gully erosion density in the black soil region has increased markedly. The most significant changes in gully erosion densities occur when cultivated land and other / degraded land were transformed. Cultivated land has somewhat lower gully density values than other/degraded but greater gully area, which demonstrated that gully erosion density on cultivated land was greatest. Cultivated land, especially on slopes is the main land use type where erosion gullies occur; (ii) The conversion of forest land and wetland into cultivated land is the dominant cause of gully erosion occur during 1965-2015. Due to the implementation of the project of "Grain for Green, " the area of forest land increased slightly since 2005, but the increasing rate of gully erosion density has not reduced even after cultivated land was converted to forest land. This seems to indicate that the conversion of cultivated land to forest has not attenuated gully erosion effectively in the recent past; (iii) Cultivating land on steep slopes may rapidly result in gully erosion, and it is more difficult to control in such situations. Erosion gullies form more easily on cultivated land and, in the case of the Kedong County, slopes greater than ∼4 • appear to be especially prone to erosion. Gully erosion control efforts should therefore initially target such slopes in the black soil region.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/supplementary material, further inquiries can be directed to the corresponding author/s.

AUTHOR CONTRIBUTIONS
SZ and TL performed field investigation. TL and ML organized the data and performed the analysis. ML wrote the original draft of the manuscript. LZ and MM reviewed and edited the manuscript. TL and WZ performed the visualization. All authors contributed to manuscript revision, read, and approved the submitted version.