# Simulating the vertical transition of soil textural layers in north-western China with a Markov chain model.

IntroductionIn the inland river basin of north-western China, oases of different sizes and shapes are commonly interspersed with widespread deserts (Cheng et al. 1999). Oases are formed naturally in the inland fiver deltas or are established on alluvial--diluvial plains and the edges of alluvial-diluvial fans (Pan 2001; Zhang et al. 2003). The textural layers of the soil profiles differ due to differing environments of deposition. Desertification in the oasis--desert ecotone (Han and Meng 1999) and the subsequent ecological restoration of desertified areas (Su et al. 2007) both influence the mechanical composition of soils. In such regions, layers of sand and clay are usually interspersed, and some profiles have impermeable layers, which significantly influence the growth of roots and the movement of water and solutes. Understanding the spatial heterogeneity and distribution of categorical soil variables (such as soil textural layers) is crucial for soil management and environmental research (Kite and Kouwen 1992; Zhu and Mackay 2001; Bouma et al. 2002). Detailed descriptions of soil profiles are usually only obtained for a small section of larger regions, because extensive sampling of profiles is impractical. Appropriate models must therefore be developed to accurately describe the vertical variability of textures in soil profiles. These models would be of great value for determining the preferential pathways of water and solute transport (Koltermann and Gorelick 1996; Weigand et al. 2001).

Geostatistics have been used to characterise the spatial heterogeneity of various soil properties. Effective methods for simulating fields of cross-correlated multinomial classes are limited, because conventional geostatistics cannot easily incorporate interclass correlations among multinomial classes (Atkinson 2001; Deutsch 2006; Li and Zhang 2006) or delineate their nonlinear spatial distributions using linear estimators (Bogaert 2002). For example, applications (Goovaerts 1996; Kyriakidis and Dungan 2001; Miller and Franklin 2002) of the indicator kriging method, introduced by Journel (1983), ignored or post-processed single realisations for the interclass correlations, which led to poor reproduction of cross-variograms (Goovaerts 1996; D'Or and Bogaert 2004).

Major methods for incorporating interclass relationships into conditional simulations have recently emerged. Markov chain theory, a theory of stochastic process, has been introduced into the geosciences. It describes how likely one state is to change to another state (Shi 1992) and can be regarded as a discrete state sequence. The probability of transition from one state to another depends on the current state; previous states are irrelevant. By using the transition probability matrix (TPM), whose elements consist of probabilities of transition from one state to another, Markov chain theory can simulate multinomial classes with the incorporation of interclass relationships. This theory has been applied to predict vegetation succession (Usher 1981; Balzter 2000; Logofet and Lesnaya 2000; Benabdellah et al. 2003; Aurbacher and Dabbert 2011; Chuang et al. 2011), canopy reflectance (Kuusk 1995), and changes in land use (Aaviksoo 1995; Zhang et al. 2011).

In soil science, categorical variables such as soil types and textural layers normally exhibit strong interdependences between multinomial classes. These interdependences include complex cross-correlations, tendencies of juxtaposition, and directional asymmetry in spatial occurrence (Li and Zhang 2006). Incorporating the interdependences is crucial to effectively capture the complex patterns of soil classes. Conditional simulation methods based on Markov chains can better reflect the features of multinomial classes (Zhang and Li 2005). A second typical feature of soil variables is that the number of classes is usually large. For example, alluvial plains or watersheds can have dozens of soil types. Markov chain methods use a TPM, with cross correlations, to accommodate this large number of classes. The one-dimensional (1-D) Markov chain theory has helped to describe and model the spatial order of parcels of different soil classes (Burgess and Webster 1984a, 1984b) and the vertical changes of soil-profile textural layers in alluvial plains (Li et al. 1997, 1999). For 2-D applications, Wu et al. (2004) developed a Markov chain Monte Carlo methodology using a neighbourhood and scanning scheme to model the 2-D spatial structure of soil. Li and Zhang (2006) extended a 2-D triplex Markov chain approach, which was based on line data (Li et al. 2004), into a generalised Markov chain model to enable a conditional and multi-dimensional simulation of categorical variables from grid-point samples. Li (2007) subsequently proposed the Markov chain random-field (MCRF) theory. Previous research, however, has been mostly focussed on alluvial or farmland soils, where landscape is uniform.

In the central region of the Heihe River system, soils, landscapes, and land uses are non-uniform. Deserts are located mainly in the north of the region. Wind erosion during desertification has sorted the soil materials, removing the fine fractions and leaving coarse-textured sand (Lobe et al. 2001; Suet al. 2004). The centre of the region contains farmland and wetland. The soil-profile textures differ significantly in this region, both vertically and horizontally. For such a complex pattern of landscapes and structures of soil layers, the vertical changes of the soil-profile texture have received little attention. The textures in soil profiles to a depth of 300 cm were thus investigated at 100 sites in a rectangular area of 100[km.sup.2]. A Markov chain--log-normal distribution (MC-LN) model was constructed to simulate the textural profiles. The objectives of this study were: (1) to investigate the vertical transition of soil textural layers by direct measurement and by simulation with a 1-D Markov chain model, and (2) to evaluate the applicability of the 1-D Markov chain model for quantitatively describing the soil textural profiles in the study area.

Materials and methods

Study area

The study was conducted in the vicinity of the Linze Inland River Basin Research Station (39[degrees]21'N, 100[degrees]07'E), Chinese Ecosystem Research Network, in Linze County, Gansu Province, China. The area has a continental arid climate. Mean annual precipitation is 120mm; ~60% of the precipitation occurs from July to September, and only 3% occurs in winter. Mean annual pan-evaporation is 2390mm, and the mean annual temperature is 7.6[degrees]C, with a minimum of -27[degrees]C and a maximum of 39.1[degrees]C. The soil is typically characterised as grey-brown desert soil (Liu et al. 2010). The farms of the oases are irrigated with water from the Heihe River and with groundwater. The aquic soil of the oases and the irrigated soil of the deserts have been formed after long-term cultivation. Saline-alkali soils and aeolian soils are extensively distributed in the area. The desert vegetations consist of Haloxylon ammodendron (C. A. Mey.) Bunge, Calligonum mongolicum Turcz., Tamarix chinensis Lour., Nitraria sphaerocarpa Maxim, and Reaumuria soongorica (Pall.) Maxim. The main crops of the oases are spring maize (Zea mays L.), spring wheat (Triticum aestivum L.), tomato (Solanum lycopersicum), and sugar beet (Beta vulgaris).

Soil sampling design

The study area is ~100[km.sup.2] (100[degrees]05'32"-100[degrees]10'01"E, 39[degrees]12'35"-39[degrees]23'28"N), with a length of 20 km from north to south and a width of 5 km from east to west, and includes three towns and nine villages. The northern part of the study area includes the southern margin of the Badain Jaran Desert, which is characterised by dense, moving and denuded residual dunes. The Heihe River flows from east to west across the centre of the study area. Landscapes are composed of farms (in oases), deserts, and wetland. A regular grid of 126 sampling sites of size 1 km by 1 km was originally designed, but the presence of the fiver channel, reservoirs, an elevated watertable, roads, and buildings forced us to abandon some sites. One hundred sites were ultimately sampled to a depth of 300 cm (Fig. 1). Sample locations were geo-referenced using a Garmin GPS 60 (Garmin Ltd, Olathe, KS) with an accuracy of 3 m.

Samples were collected using a hand auger 5 cm in diameter. The total sampling depth of 300cm was divided into two subintervals: 0-100 and 100-300cm. Samples in the 0-100cm subinterval were collected every 10cm and samples in the 100-300cm subinterval were collected every 20 cm by pooling two 10-cm subsamples. All 2000 samples were air-dried, sealed in air-tight bags, and taken to the laboratory. After passing the samples through a 2-mm, stainless-steel sieve, the soil-particle composition was measured using a Malvern Laser particle-size analyser (Malvern Instruments Ltd, Malvern, UK) (Liu et al. 2005). The measured soil textures covered nine classes of the US Department of Agriculture's Soil Taxonomy (Shirazi and Boersma 1984), i.e. sand, loamy sand, sandy loam, silt loam, loam, clay loam, silty clay loam, silty clay, and clay, respectively. The textural type, thickness, and sequence of the different layers were determined.

Simulation model

Alluvial soils, distributed on the river banks, are usually composed of several profile textural layers (Fig. 1). In the vertical direction, soil profiles are space-continuous, and textural types are discrete. The vertical change of textural layers is assumed to be a 1-D Markov chain that can be described by a TPM. The calculation of the TPM in the present study is based on the actual change of textural types, because transitions among different textural types are easy to observe and describe. The layer thickness of each textural type is described by a proper probability density function.

The distribution of the raw dataset of textural layer thicknesses was skewed, so the data were logarithmically (ln) transformed and retested for normality with a one-sample Kolmogorov-Smirnov (K-S) test.

The probability density function of a lognormal distribution can be expressed as:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (1)

where x is a value of textural-layer thickness and [mu] and [sigma] are the mean and standard deviation of the lognormal distribution, respectively. For the layer thicknesses of textural types, the characteristic parameters [mu] and [sigma] can be described as:

[mu] = {[[mu].sub.i]} = ([[mu].sub.i], [[mu].sub.2], ..., [[mu].sub.m])

[sigma] = {[[sigma].sub.i]} = ([[sigma].sub.1], [[sigma].sub.2], ..., [[sigma].sub.m]), (i = 1, 2, ..., m) (2)

where m is the total number of the textural types (m = 9); i = 1,2, ..., 9 represent sand, loamy sand, sandy loam, silt loam, loam, clay loam, silty clay loam, silty clay, and clay, sequentially, and [[mu].sub.i] and [[sigma].sub.i] are the logarithmic mean and standard deviation of the thickness of the ith textural type, respectively.

The type, sequence, and thickness of the soil-profile textural layers in the study area were recorded. Transitions from top to bottom between two fixed and adjoined textural layers in the entire area were summed, and the transition frequency matrix (TFM) N was:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (3)

where m is the total number of the textural types and [n.sub.ij] is the total transition frequency from the ith textural type into the jth type. For two neighbouring profile layers, the jth type is located below the ith type; [n.sub.i] is the total transition frequency outward from the ith type ([n.sub.i] = [[summation].sup.m.sub.j=1] [n.sub.ij]).

Based on the transition frequency [n.sub.ij], the transition probability [p.sub.ij] can be obtained by dividing [n.sub.ij] by [n.sub.i], so the TPM (P) is calculated as:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (4)

In the TPM, the diagonal elements [p.sub.ii] (auto-transitions, [p.sub.ii] = 0) represent autocorrelations of individual classes, whereas the off-diagonal elements [p.sub.ij](i[not equal to]j) (cross-transitions) represent cross-correlations between different classes (Li et al. 2005).

The TPM combined with the initial probability distribution (IPD) can completely describe the transition characteristics of different textural types. The initial probability equals the ratio of the area of a certain surface textural type in the entire study area (Li et al. 1999):

A = {[A.sub.i]} = ([A.sub.1], [A.sub.2], ..., [A.sub.m]), (i = 1,2 ... m) (5)

where [A.sub.i] is the initial probability of the ith textural type and m is the total number of the textural types.

Using a Monte Carlo simulation, the accumulative IPD (CA) and the accumulative TPM (CP) are required:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (6)

where [CA.sub.i] is the accumulative initial probability from the first textural type (sand) to the ith type ([CA.sub.i] = [A.sub.1] + [A.sub.2] + ... + [A.sub.i]), [CP.sub.ij] is the accumulative probability for the ith textural type transitioning from the first type (sand) to the jth type ([CP.sub.ij] = [P.sub.i1] + [P.sub.i2] + ... + [P.sub.ij], (i, j = l, 2 ...., m)), and m is the total number of the textural types.

Test of the Markov chain characteristic

Transitions between neighbouring textural layers may be correlated and assumed to be a 1-D Markov chain. The existence of correlations and the Markov characteristic, however, must be verified. The null hypothesis is that changes of textural layers are independent of each other. Otherwise, a Markov chain is assumed. Anderson and Goodman (1957) proposed the following statistical variable to test:

-2ln[lambda] = 2 [m.summation over (i,j=1)][n.sub.ij] In ([p.sub.ij]/[p.sub.j]) (7)

where [lambda], is the statistical variable, [p.sub.j] is the edge probability of the jth column in the TPM ([p.sub.j] = [[sigma].sup.m.sub.i=1] [n.sub.ij]/[[sigma].sup.m.sub.ij=1][n.sub.ij]) and [n.sub.ij],[p.sub.ij], and m are as above.

If the null hypothesis is true, -2ln[lambda] approximates a [chi square] distribution with a degree of freedom (d.f.)= m(m - 1).

Test of the stability of vertical textural-layer changes

In our study, [p.sub.ij] is the probability that textural layer i at depth x - 1 transitions vertically to type j at depth x. If the Markov chain is stationary for the entire study area, [p.sub.ij]. is irrelevant to the position x, ([p.sub.ij](x)=[p.sub.ij]). The entire region or the total depth is divided into several subintervals, and the corresponding TPMs are calculated. If the TPMs of the subintervals agree well with that of the whole, the TPM of the whole can reflect the vertical change of the profile textural layers. Otherwise, [p.sub.ij] is a function of position x, and it is inappropriate to reflect the textural-layer transitions using the TPM of the whole. Anderson and Goodman (1957) presented the following criterion to test the stability:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (8)

where [lambda] is the statistical variable, [p.sub.ij.sup.(t)] is the (i,j)th element in the TPM of the tth subinterval, [p.sub.ij] is the (i,j)th element in the TPM of the whole, [n.sub.ij.sup.(t)] is the (i, j)th element in the TFM of the tth subinterval, and T and m are the total number of the subintervals and the textural types, respectively.

When the null hypothesis, i.e. a stationary Markov chain, is true, -2ln[lambda] accords with a [chi square] distribution with d.f. = (T - 1) [m(m - 1)]. Otherwise, if the test statistic exceeds the corresponding [chi square] value, the null hypothesis is rejected, and the test for the order of the chain must be performed (Balzter 2000; Weigand et al. 2001).

Simulation

When the aforementioned verifications were satisfactory, the MC-LN model constructed based on the Monte Carlo simulation (Li et al. 1999), using Matlab 7.0 software (www.mathworks.com), was applied to simulate the textural profiles (Fig. 2).

Criteria to test the results of simulation

The TPM was calculated based on the simulated soil profiles, and the relative error (RE) was used as a criterion (Hu et al. 2012) to evaluate the deviations between the measured and the simulated TPM:

RE = [absolute value of ([p'.sub.ij] - [p.sub.ij])]/[p.sub.ij] x 100% (10)

where [p'.sub.ij] and [p.sub.ij] are the (i,j)th element in the simulated and the measured TPM, respectively.

Results and discussion

Probability distribution of textural-layer thickness

The descriptive statistics of the textural-layer thicknesses are summarised in Table 1. The mean layer thickness of each textural type differed, with relatively thick layers of sand and silty clay loam, and relatively thin layers of sandy loam, loam, and clay. The extensive distribution of deserts in the northern and central parts of the study area produced distinctly thicker (means) layers of sand relative to the other textural types. The occurrence of layers of silt loam and clay was relatively rare, which consequently influenced the mean thicknesses. The coefficients of variation (CV) ranged from 36% to 87%, indicating moderate variation of textural-layer thicknesses (Nielsen and Bouma 1985), among which the variation of silt loam was relatively weak (CV = 36%) while that of silty clay loam was relatively strong (CV = 87%).

The K-S test indicated that the layer thicknesses of each textural type were log-normally distributed (P>0.05). This conclusion is consistent with the report of Li et al. (1997). The characteristic parameters [mu] and [sigma] fit by the log-normal distribution were:

[mu] = (4.58, 3.73, 3.60, 4.10, 3.63, 3.70, 4.48, 3.93, 3.59)

[sigma] = (0.82, 0.76, 0.76, 0.36, 0.73, 0.78, 0.76, 0.84, 0.76)

Transition characteristics of different textural layers

Due to difficulties in surveying the areas of surface (0-10 cm) textural types, the elements in the IPD were expressed as the occurrence probabilities of different surface textural types (Li et al. 1999) in the 100 profiles as:

A = (0.45, 0.02, 0.06, 0.01, 0.15, 0.17, 0.14, 0.00, 0.00)

The accumulative IPD was:

CA = (0.45, 0.47, 0.53, 0.54, 0.69, 0.86, 1, 1, 1)

In the IPD, [A.sub.8] = [A.sub.9] = 0, while in the accumulative IPD, [CA.sub.7] = [CA.sub.8] = [CA.sub.9] = 1. Excluding the layers of silty clay and clay, the other seven types of textural layer all occurred in the surface soil with different probabilities. Sand was the most common textural type in the surface soil.

The TFM (N), TPM (P), and the accumulative TPM (CP) of different textural types were calculated:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Test of the Markov chain characteristic

A test of the Markov-chain characteristic indicated that [[chi square].sub.0.05] (70) = 90.53, which was much lower than -2ln[lambda] = 274.12, with d.f. = m(m - 1) = 72. The alterative hypothesis that transitions between neighbouring textural layers had some correlations was thus accepted, and the vertical transition of different textural layers was regarded as a 1-D Markov chain.

Test of the stability of the vertical transition of textural layers

Environments and land uses conducive to sedimentation are vital for the formation of soil textural layers. The types, sequences, and thicknesses of profile textural layers differ as the depositing environment changes. The occurrence of a certain textural layer or of a fixed combination of textural layers more than twice in a prolife are rare. For example, profiles in the study area contained one to seven textural types, and layers of sand dominated the 100 profiles. Profiles were composed solely of layers of sand at 11 sites, and four profiles were composed solely of layers of silty clay loam. Testing the stability of the vertical transition of soil textural layers was thus necessary.

The corresponding TPMs of the 0-100 and 100-300cm subintervals were then calculated (Tables 2 and 3). More layers occurring below a specified textural layer in the 10-300cm subinterval than in the 0-100cm subinterval indicated the more obvious layering structure for the 100-300cm subinterval. On the other hand, fewer layers occurred beneath layers of loamy sand and silt loam in both subintervals. Layers of silty clay loam occurred beneath layers of clay loam and silty clay, and layers of silty clay occurred beneath layers of silty clay loam, with relatively high probabilities in both subintervals. The different sequences of soil textural layers in the 0-100 and 100-300cm subintervals may have been caused by climate change, geological movements, or anthropogenic activity during the processes of soil formation. The migration of the Heihe River over time may also have significantly influenced the soil texture of the beaches. The statistical analysis indicated that -2ln[lambda] = 83.31, with d.f. = (2 - 1) x 9 x (9 - 1) = 72, was slightly less than [[chi square].sub.0.05](70)=90.53. The null hypothesis of no significant difference in the elements of the TPMs between the entire profile and the subintervals was accepted, so the vertical transition of soil textural layers could be considered as stable.

Because the soils and landscapes on the northern and southern banks of the Heihe River differed, the entire study area was divided into northern and southern subregions bounded by the Heihe River. Differences and similarities were observed between the TPMs of the two subregions (Tables 4 and 5). The main difference between the two subregions was that textural layers beneath layers of silt loam and clay were more diverse in the southern subregion than in the northern subregion. Additionally, transitions between fine-textured soil layers were more common in the southern subregion. These differences were likely due to different sedimentary environments and conditions of soil genesis. The extensively distributed wetland in the southern subregion originated from snowmelt from the Qilian Mountains. Springs here overflow in lowlands, and seasonal or perennial stagnant water affects the formation of soil. Soil types are meadow Solonchak (FAO 1998) and bog with high contents of fine particles. The northern subregion, however, includes the southern margin of the Badain Jaran Desert. Careless use of water and soil has degraded the oases and desertified the land around them, which is now characterised by coarse-textured soils. The statistical analysis indicated that -2ln[lambda] = 62.20, with d.f. =(2 - 1) x 9 x (9 - 1) = 72, was less than [[chi square].sub.0.05](70) = 90.53. The hypothesis that no significant difference existed between the TPMs of the subregions and entire area was accepted, so the vertical transition of soil textural layers in the entire study area was stable.

Based on the relatively high element values in the TPMs of the subintervals, subregions, and entire study area, the main combinations of upper to lower textural layers in the soil profiles were loamy sand and sand (or sandy loam), sandy loam and sand (or loamy sand and sand), loam and clay loam, clay loam (or silty clay) and silty clay loam, and silty clay loam and silty clay.

Analysis of the simulated layer thicknesses of textural types

In the MC-LN model based on Monte Carlo simulation, textural layers were produced using a sequence of random numbers, and layer thicknesses were produced using log-normally distributed random numbers. For each textural type, the simulated mean layer thickness agreed well with the measured value (Table 6). The simulated mean layer thickness of clay loam was the closest to the measured value (RE = 8%). The relative deviations between the measured and the simulated mean layer thickness of silty clay loam (RE=24%), sand (RE=23%), and silty clay (RE = 22%) were larger than those of the other textural types.

The differences between the measured and the simulated mean layer thicknesses could be ascribed to two factors. Errors are introduced when calculating the thicknesses of the simulated textural layers. During sampling, soil texture within the interval of sampling depth (10cm for the 0-100cm subinterval and 20cm for the 100-300cm subinterval) was considered homogeneous, and thinner layers (<10 or 20cm) were ignored due to their rare occurrence. The simulated layer thicknesses, though, were produced based on the log-normally distributed random numbers in the Monte Carlo simulation, and the simulated textural layers included thinner (<10 or 20cm) layers. The simulated mean layer thickness of each textural type may thus be slightly smaller than the measured values. The other factor is bottom-layer cutoff (Li et al. 1999). When the bottom layer was so thick that the actual profile depth exceeded 300cm, the bottom layer was 'cut off', especially when the bottom layer was immediately transitioned from other textural types. This cutoff produced a higher probability of occurrence of thin layers in the bottom, and the simulated mean layer thicknesses tended to be slightly smaller than the measured values. The frequencies of occurrence of layers of silty clay loam, sand, and silty clay in the bottom were 20, 38, and 9, respectively, which were relatively higher than those of other types. For example, layers of clay loam only occurred in the bottom two times, while the simulated mean layer thickness was the closet to the measured value. This evidenced the influence of bottom layer cutoff on the simulated mean layer thicknesses of the textural types.

Analysis of the simulated vertical change of textural layers

Based on the simulated soil profiles, the simulated TPM ([P.sub.s]) was given as:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

The occurrences of layers of loamy sand, silt loam, and clay in the study area were rare (Table 1), so their occurrences and sequences in profiles were much more random and heterogeneous. Consequently, the transitions from other textural types into these three types were less frequent. None of these three layers occurred mutually beneath each other. The simulated TPM approximately equalled the measured TPM (Table 7). The largest relative deviation (RE = 12.6%) occurred when layers of clay loam transitioned into layers of clay. The second (10.4%) and the third (10.2%) largest relative deviations occurred when layers of silt loam transitioned into sandy loam and when layers of sandy loam transitioned into layers of clay loam, respectively. Li et al. (1999) pointed out that the random numbers created by Monte Carlo simulation were 'false random numbers'. For a given random number seed, the sequence of random numbers produced was unique. Deviations thus existed between distributions of random number sequences produced from different seeds. We, however, used Matlab 7.0 software to run the simulations. Random numbers produced by Matlab 7.0 are nearly 'true random numbers', and the simulative precision of the MC-LN model is improved. On the other hand, the transition of different textural layers was verified as stationary, i.e. the probabilities of transition did not change with position. The probabilities of transition, though, might not remain constant, but may increase or decrease during the transition from one state to another (Usher 1981; Li et al. 1997). The TPM is thus only approximately stationary. Differences exist between the TPMs of the subintervals and the profiles and between the TPMs of the subregions and the entire study area. The simulated results of the stationary Markov chain may thus deviate from the measured data. Because the profile depth and the number of textural layers are limited, the occurrence of the same textural layer in a profile is rarely repeated, and the application of the stationary Markov chain is convenient and is easy to implement, so the stationary Markov chain is generally acceptable. The stationary MC-LN model in our study quantitatively simulated the vertical transition of soil-profile textures and layer thicknesses.

Conclusions

A quantitative description of the vertical transition of soil textural layers is important for quantifying the transport of water and solutes. After investigating the vertical variability of soil-profile textures to a depth of 300cm at 100 sites, an MC-LN model was successfully constructed and used to simulate the soil textural profiles. We concluded that the vertical transition of adjoined textural layers in soil profiles in the study area could be regarded as a 1-D Markov chain. The vertical transition of different textural types in the entire study area was stationary. The MC-LN model accurately simulated the soil textural profiles. The main combinations of upper to lower textural layers in the soil profiles were loamy sand and sand (or sandy loam), sandy loam and sand (or loamy sand and loam), loam and clay loam, clay loam (or silty clay) and silty clay loam, and silty clay loam and silty clay.

In the central region of the Heihe River system, interchanges between soil water and groundwater in farmland and wetland frequently occur owing to the hydraulic connections. Severe evaporation, irrigation of farmland by groundwater rich in salts, and seasonal stagnant water in the wetland result in salinisation. Integrating the vertical heterogeneity of soil texture modelled by Markov chains with models of water and solute transport will be another avenue for gaining insight into the characteristics of soil water consumption (Ye and Khaleel 2008; Sivakumar et al. 2005; Sun and Li 2010). This could be a new tactic for learning how to improve water-use efficiency and how to relieve salinisation in arid regions in north-western China.

http://dx.doi.org/10.1071/SR12332

Acknowledgements

This study was financially supported by the National Natural Science Foundation of China (Project number: 91025018). Authors express our sincere appreciation to editors and two anonymous reviewers for their useful comments and suggestions. Mr Robert Horton is greatly appreciated for checking the English of this manuscript. Thanks go to staff of the State Key Laboratory of Soil Erosion and Dryland Farming on the Loess Plateau, Institute of Soil and Water Conservation, Chinese Academy of Sciences.

Received 14 November 2012, accepted 27 May 2013, published online 18 June 2013

References

Aaviksoo K (1995) Simulating vegetation dynamics and land use in a mire landscape using a Markov model. Landscape and Urban Planning 31, 129-142. doi:10.1016/0169-2046(94)01045-A

Anderson TW, Goodman LA (1957) Statistical inference about Markov chains. Annals of Mathematical Statistics 28, 89-110. doi:10.1214/aoms/l177707039

Atkinson PM (2001) Geographical information science: GeoComputation and nonstationarity. Progress in Physical Geography 25, 111-122.

Aurbacher J, Dabbert S (2011) Generating crop sequences in land-use models using maximum entropy and Markov chains. Agricultural Systems 104, 470-479. doi:10.1016/j.agsy.2011.03.004

Balzter H (2000) Markov chain models for vegetation dynamics. Ecological Modelling 126, 139-154. doi:10.1016/S0304-3800(00)00262-3

Benabdellah B, Albrecht KF, Pomaz VL, Denisenko EA, Logofet DO (2003) Markov chain models for forest successions in Erzgebirge, Germany. Ecological Modelling 159, 145-160. doi:10.1016/S0304-3800(02)00285-5

Bogaert P (2002) Spatial prediction of categorical variables: the Bayesian maximum entropy approach. Stochastic Environmental Research and Risk Assessment 16, 425-448. doi:10.1007/s00477-002-0114-4

Bouma J, van Alphen BJ, Stoorvogel JJ (2002) Fine tuning water quality regulations in agriculture to soil differences. Environmental Science & Policy 5, 113-120. doi:10.1016/S1462-9011(02)00030-8

Burgess TM, Webster R (1984a) Optimal sampling strategies for mapping soil types. 1. Distribution of boundary spacings. Journal of Soil Science 35, 641-4554. doi:10.111 l/j.1365-2389.1984.tb00621.x

Burgess TM, Webster R (1984b) Optimal sampling strategies for mapping soil types, II. Risk functions and sampling intervals. Journal of Soil Science 35, 655-665. doi:10.111l/j.1365-2389.1984.tb00622.x

Cheng GD, Xiao DN, Wang GX (1999) On the characteristics and building of landscape ecology in arid area. Advances in Earth Science 14, 11-15. [in Chinese with English abstract]

Chuang CW, Lin CY, Chien CH, Chou WC (2011) Application of Markov-chain model for vegetation restoration assessment at landslide areas caused by a catastrophic earthquake in Central Taiwan. Ecological Modelling 222, 835-845. doi:10.1016/j.ecolmodel.2010.11.007

D'Or D, Bogaert P (2004) Spatial prediction of categorical variables with the Bayesian maximum entropy approach: the Ooypolder case study. European Journal of Soil Science 55, 763-775. doi:10.1111/j.1365-2389.2004.00628.x

Deutsch CV (2006) A sequential indicator simulation program for categorical variables with point and block data: BlockSIS. Computers & Geosciences 32, 1669-1681. doi:10.1016/j.cageo.2006.03.005

FAO (1998) World reference base for soil resources. In 'Key to the reference soil groups'. World Soil Resources Report 84. (FAO: Rome)

Goovaerts P (1996) Stochastic simulation of categorical variables using a classification algorithm and simulated annealing. Mathematical Geology 28, 909-921. doi:l0.1007/BF02066008

Han DL, Meng XY (1999) Recent progress of research on oasis in China. Chinese Geographical Science 9, 199-205. doi:10.1007/s11769-999-0044-x

Hu W, Tallon LK, Si BC (2012) Elevation of time stability indices for soil water storage upscaling. Journal of Hydrology 475, 229-241. doi:10.1016/j.jhydrol.2012.09.050

Journel AG (1983) Nonparametric estimation of spatial distributions. Journal of the International Association for Mathematical Geology 15, 445-468. doi:10.1007/BF01031292

Kite GW, Kouwen N (1992) Watershed modeling using land classifications. Water Resources Research 28, 3193-3200. doi:10.1029/92WR01819

Koltermann CE, Gorelick SM (1996) Heterogeneity in sedimentary deposits: a review of structure-imitating, process-imitating, and descriptive approaches. Water Resources Research 32, 2617-2658. doi:10.1029/96WR00025

Kuusk A (1995) A Markov chain model of canopy reflectance. Agricultural and Forest Meteorology 76, 221-236. doi:10.1016/0168-1923(94)02216-7

Kyriakidis PC, Dungan JL (2001) A geostatistical approach for mapping thematic classification accuracy and evaluating the impact of inaccurate spatial data on ecological model predictions. Environmental and Ecological Statistics 8, 311-330. doi:10.1023/A:1012778302005

Li WD (2007) Markov chain random fields for estimation of categorical variables. Mathematical Geology 39, 321-335. doi:10.1007/s11004-007-9081-0

Li WD, Zhang CR (2006) A generalized Markov chain approach for conditional simulation of categorical variables from grid samples. Transactions in GIS 10, 651-669. doi:10.1111/j.1467-9671.2006.01017.x

Li WD, Li BG, Shi YC, Tang DY (1997) Application of the Markov chain theory to describe spatial distribution of textural layers. Soil Science 162, 672-683. doi: 10.1097/00010694-199709000-00009

Li WD, Li BG, Shi YC (1999) Markov-chain simulation of soil textural profiles. Geoderma 92, 37-53. doi: 10.1016/S0016-7061(99)00024-5

Li WD, Zhang CR, Burt JE, Zhu AX, Feyan J (2004) Two-dimensional Markov chain simulation of soil type spatial distribution. Soil Science Society of America Journal 68, 1479-1490. doi:10.2136/sssaj2004.1479

Li WD, Zhang CR, Burt JE, Zhu AX (2005) A Markov chain-based probability vector approach for modeling spatial uncertainties of soil classes. Soil Science Society of America Journal 69, 1931-1942.

Liu YP, Tong J, Li XN (2005) Analysing the silt particles with the Malvern Mastersizer 2000. Water Conservancy Science and Technology and Economy 11, 329-331. [in Chinese with English abstract]

Liu B, Zhao WZ, Chang XX, Li SB, Zhang ZH, Du MW (2010) Water requirements and stability of oasis ecosystem in arid region, China. Environmental Earth Science 59, 1235-1244. doi:10.1007/s12665-009-0112-7

Lobe I, Amelung W, Du Preez CC (2001) Losses of carbon and nitrogen with prolonged arable cropping from sandy soils of the South Africa Highveld. European Journal of Soil Science 52, 93-101. doi:10.1046/j.1365-2389.2001.t01-1-00362.x

Logofet DO, Lesnaya EV (2000) The mathematics of Markov models: what Markov chains can really predict in forest successions. Ecological Modelling 126, 285-298. doi:10.1016/S0304-3800(00)00269-6

Miller J, Franklin J (2002) Modeling the distribution of four vegetation alliances using generalized linear models and classification trees with spatial dependence. Ecological Modelling 157, 227-247. doi:10.1016/S0304-3800(02)00196-5

Nielsen DR, Bouma J (1985) Soil spatial variability. In 'Proceedings of Workshop of the ISSS and the SSSA'. 30 Nov.-1 Dec. 1984, Las Vegas, NV. (Pudoc: Wageningen, the Netherlands)

Pan XL (2001) A preliminary study on the stability of oasis ecosystem in arid area. Quaternary Sciences 21, 345-351. [in Chinese with English abstract]

Shi RJ (1992) 'The basic theory of Markov chain and its application.' (Xi'an Electronic and Technical University Publishing: Xi'an, China) [in Chinese].

Shirazi MA, Boersma L (1984) A unifying quantitative analysis of soil texture. Soil Science Society of America Journal 48, 142-147. doi:10.2136/sssaj1984.03615995004800010026x

Sivakumar B, Hatter T, Zhang H (2005) A fractal investigation of solute travel time in a heterogeneous aquifer: transition probability/Markov chain representation. Ecological Modelling 182, 355-370. doi:10.1016/j.ecolmodel.2004.04.010

Su YZ, Zhao HL, Zhao WZ, Zhang TH (2004) Fractal features of soil particle size distribution and the implication for indicating desertification. Geoderma 122, 43-49. doi:10.1016/j.geoderma.2003.12.003

Su YZ, Zhao WZ, Su PX, Zhang ZH, Wang T, Ram R (2007) Ecological effects of desertification control and desertified land reclamation in an oasis-desert ecotone in an arid region: A case study in Hexi Corridor, northwest China. Ecological Engineering 29, 117-124. doi:10.1016/j.ecoleng.2005.10.015

Sun DZ, Li XQ (2010) Application of Markov Chain model on environmental fate of phenanthrene in soil and groundwater. Procedia Environmental Science 2, 814-823. doi:10.1016/j.proenv.2010.10.092

Usher MB (1981) Modelling ecological succession, with particular reference to Markovian models. Vegetation 46-47, 11-18. doi:10.1007/BF00118380

Weigand H, Totsche KU, Huwe B, Kogel-Knabner I (2001) PAH mobility in contaminated industrial soils: a Markov chain approach to the spatial variability of soil properties and PAH levels. Geoderma 102, 371-389. doi:10.1016/S0016-7061(01)00043-X

Wu KJ, Nunan N, Crawford JW, Young IM, Ritz K (2004) An efficient Markov chain model for the simulation of heterogeneous soil structure. Soil Science Society of America Journal 68, 346-351. doi:10.2136/sssaj2004.0346

Ye M, Khaleel R (2008) A Markov chain model for characterizing medium heterogeneity and sediment layering structure. Water Resources Research 44, W09427. doi:10.1029/2008WR006924

Zhang CR, Li WD (2005) Markov chain modeling of multinomial land-cover classes. GISeience & Remote Sensing 42, 118. doi:10.2747/1548-1603.42.1.1

Zhang H, Wu JW, Zheng QH, Yu YJ (2003) A preliminary study of oasis evolution in the Tarim Basin, Xinjiang, China. Journal of Arid Environments 55, 545-553. doi:10.1016/S0140-1963(02)00283-5

Zhang RQ, Tang CJ, Ma SH, Yuan H, Gao LL, Fan WY (2011) Using Markov chains to analyze changes in wetland trends in arid Yinchuan Plain, China. Mathematical and Computer Modelling 54, 924-930. doi:10.1016/j.mcm.2010.11.017

Zhu AX, Mackay DS (2001) Effects of spatial detail of soil information on watershed modeling. Journal of Hydrology 248, 54-77. doi:10.1016/S0022-1694(01)00390-0

Danfeng Li (A,C) and Ming'an Shao (B,D)

(A) State Key Laboratory of Soil Erosion and Dryland Farming on the Loess Plateau, Institute of Soil and Water Conservation, Chinese Academy of Sciences and Ministry of Water Resources, Yangling, 712100 Shaanxi, P.R. China.

(B) Key Laboratory of Ecosystem Network Observation and Modeling, Institute of Geographic Sciences and Natural Resources Research, Chinese Academy of Sciences, Beijing 100101, P.R. China.

(C) University of Chinese Academy of Sciences, Beijing 100049, P.R. China.

(D) Corresponding author. Email: shaoma@igsnrr.ac.cn

Table 1. Descriptive statistics of the layer thicknesses of the textural types in soil profiles of the study area Frequency is the total occurrence number of each textural type with thickness ranging from 10 to 300 cm at the 100 sites. s.d., Standard deviation; CV, coefficient of variation Statistic Sand Loamy Sandy Silt Loam Clay sand loam loam loam Frequency 83 26 44 12 45 53 Mean (cm) 131 55 49 64 48 53 Max. (cm) 300 200 230 100 140 180 Min. (cm) 20 10 10 40 10 10 s.d. (cm) 95 43 41 23 34 39 CV (%) 72 79 84 36 70 73 Statistic Silty clay Silty Clay loam clay Frequency 58 36 13 Mean (cm) 114 72 45 Max. (cm) 300 250 100 Min. (cm) 10 10 10 s.d. (cm) 78 63 29 CV (%) 69 87 63 Table 2. Transition probability matrix of the textural types in the 0-100 cm subinterval in the study area Upper layer Lower layer Sand Loamy Sandy Silt Loam sand loam loam Sand 0 0.30 0.25 0.00 0.05 Loamy sand 0.56 0 0.44 0.00 0.00 Sandy loam 0.15 0.31 0 0.00 0.38 Silt loam 0.00 0.00 1.00 0 0.00 Loam 0.09 0.09 0.00 0.00 0 Clay loam 0.04 0.07 0.07 0.04 0.11 Silty clay loam 0.00 0.10 0.00 0.20 0.00 Silty clay 0.00 0.00 0.13 0.00 0.13 Clay 0.00 0.00 0.00 0.00 0.00 Upper layer Clay Silty clay Silty Clay loam loam clay Sand 0.00 0.10 0.20 0.10 Loamy sand 0.00 0.00 0.00 0.00 Sandy loam 0.15 0.00 0.00 0.00 Silt loam 0.00 0.00 0.00 0.00 Loam 0.64 0.09 0.09 0.00 Clay loam 0 0.48 0.15 0.04 Silty clay loam 0.30 0 0.40 0.00 Silty clay 0.13 0.50 0 0.13 Clay 0.00 0.00 1.00 0 Table 3. Transition probability matrix of the textural types in the 100-300 cm subinterval in the study area Upper layer Lower layer Sand Loamy Sandy Silt Loam sand loam loam Sand 0 0.04 0.24 0.00 0.14 Loamy sand 0.36 0 0.27 0.00 0.27 Sandy loam 0.47 0.16 0 0.05 0.11 Silt loam 0.40 0.00 0.00 0 0.60 Loam 0.11 0.00 0.33 0.17 0 Clay loam 0.19 0.10 0.10 0.05 0.05 Silty clay loam 0.15 0.07 0.04 0.07 0.11 Silty clay 0.05 0.00 0.05 0.05 0.05 Clay 0.14 0.00 0.29 0.00 0.14 Upper layer Clay Silty Clay Silty Clay loam loam clay Sand 0.14 0.19 0.05 0.19 Loamy sand 0.00 0.09 0.00 0.00 Sandy loam 0.05 0.05 0.05 0.05 Silt loam 0.00 0.00 0.00 0.00 Loam 0.22 0.11 0.06 0.00 Clay loam 0 0.33 0.19 0.00 Silty clay loam 0.07 0 0.48 0.00 Silty clay 0.05 0.32 0 0.42 Clay 0.29 0.00 0.14 0 Table 4. Transition probability matrix of the textural types in the northern subregion in the study area Upper layer Lower layer Sand Loamy Sandy Silt Loam sand loam loam Sand 0 0.23 0.31 0.00 0.08 Loamy sand 0.62 0 0.31 0.00 0.08 Sandy loam 0.33 0.27 0 0.00 0.27 Silt loam 0.00 0.00 0.00 0 1.00 Loam 0.11 0.05 0.16 0.05 0 Clay loam 0.13 0.13 0.00 0.00 0.06 Silty clay loam 0.14 0.07 0.00 0.14 0.14 Silty clay 0.17 0.00 0.00 0.00 0.17 Clay 1.00 0.00 0.00 0.00 0.00 Upper layer Clay Silty clay Silty Clay loam loam clay Sand 0.04 0.23 0.08 0.04 Loamy sand 0.00 0.00 0.00 0.00 Sandy loam 0.07 0.00 0.07 0.00 Silt loam 0.00 0.00 0.00 0.00 Loam 0.42 0.16 0.05 0.00 Clay loam 0 0.69 0.00 0.00 Silty clay loam 0.00 0 0.50 0.00 Silty clay 0.00 0.50 0 0.17 Clay 0.00 0.00 0.00 0 Table 5. Transition probability matrix of the textural types in the southern subregion in the study area Upper layer Lower layer Sand Loamy Sandy Silt Loam sand loam loam Sand 0 0.10 0.19 0.00 0.14 Loamy sand 0.22 0 0.44 0.00 0.22 Sandy loam 0.35 0.15 0 0.05 0.15 Silt loam 0.29 0.00 0.14 0 0.57 Loam 0.10 0.05 0.14 0.10 0 Clay loam 0.09 0.06 0.15 0.06 0.09 Silty clay loam 0.09 0.08 0.08 0.08 0.04 Silty clay 0.00 0.00 0.13 0.04 0.09 Clay 0.11 0.00 0.22 0.00 0.11 Upper layer Clay Silty clay Silty Clay loam loam clay Sand 0.09 0.05 0.19 0.24 Loamy sand 0.00 0.11 0.00 0.00 Sandy loam 0.20 0.05 0.00 0.05 Silt loam 0.00 0.00 0.00 0.00 Loam 0.48 0.05 0.10 0.00 Clay loam 0 0.29 0.24 0.03 Silty clay loam 0.21 0 0.42 0.00 Silty clay 0.09 0.30 0 0.35 Clay 0.22 0.00 0.33 0 Table 6. Comparisons between the measured and simulated mean layer thicknesses (cm) of the textural types RE, Relative error Sand Loamy Sandy Silt Loam sand loam loam Measured 131 55 49 64 48 Simulated 101 48 43 56 43 RE (%) 23 12 12 12 10 Clay Silty clay Silty Clay loam loam clay Measured 53 114 72 45 Simulated 49 87 56 40 RE (%) 8 24 22 11 Table 7. Per cent relative error values (%) between the measured and simulated transition probability matrix Upper layer Lower layer Sand Loamy Sandy Silt Loam sand loam loam Sand 0 0.7 2.1 0.0 3.2 Loamy sand 2.1 0 1.8 0.0 0.7 Sandy loam 1.9 0.9 0 9.0 3.3 Silt loam 0.7 0.0 10.4 0 1.8 Loam 3.2 6.6 0.1 8.8 0 Clay loam 0.3 0.6 3.2 2.2 5.7 Silty clay loam 2.4 0.4 1.2 0.5 1.8 Silty clay 1.1 0.0 7.8 2.1 2.4 Clay 5.6 0.0 5.6 0.0 0.6 Upper layer Clay Silty clay Silty Clay loam loam clay Sand 8.7 1.3 2.9 3.4 Loamy sand 0.0 8.3 0.0 0.0 Sandy loam 10.2 1.7 2.0 3.9 Silt loam 0.0 0.0 0.0 0.0 Loam 0.8 8.5 6.9 0.0 Clay loam 0 0.2 5.4 12.6 Silty clay loam 2.7 0 0.7 0.0 Silty clay 5.9 2.3 0 2.0 Clay 5.0 0.0 3.5 0

Printer friendly Cite/link Email Feedback | |

Author: | Li, Danfeng; Shao, Mingan |
---|---|

Publication: | Soil Research |

Article Type: | Report |

Geographic Code: | 9CHIN |

Date: | May 1, 2013 |

Words: | 7730 |

Previous Article: | The World Reference Base for soils (WRB) and Soil taxonomy: an appraisal of their application to the soils of the Northern Rivers of New South Wales. |

Next Article: | Preconsolidation pressure, soil water retention characteristics, and texture of Latosols in the Brazilian Cerrado. |

Topics: |