Fractal Characteristics and Prediction of Backsilting Quantity in Yangtze Estuary Deepwater Channel

Fractal interpolation has been an important method applied to engineering in recent years. It can not only be used to fit smooth curve and stationary data but also show its unique superiorities in the fatting of non-smooth curve and non-stationary data. Through analyzing such characteristic values as average value, standard deviations, skewness and kurtosis of measured backsilting quantities in the Yangtze Estuary 12.5 m Deepwater Channel during 2011–2017, the fractal interpolation method can be used to study the backsilting quantity distribution with time. According to the fractal interpolation made on the channel backsilting quantities from January 2011 to December 2017, there was a good corresponding relationship between the annual (monthly) siltation quantities and the vertical scaling factor. On this basis, a calculation formula for prediction of the backsilting quantity in the Yangtze Estuary Deepwater Channel was constructed. With the relationship between the predicted annual backsilting quantities and the vertical scaling factor, the monthly backsilting quantities can be obtained. Thus, it provides a new method for estimating the backsilting quantity of the Yangtze Estuary Deepwater Channel.


Introduction
The Yangtze River is the third largest river in the world, entering into the East China Sea from a nearby area of Shanghai, China. The Yangtze Estuary is a silty and muddy estuary with three-stage branches and four entrances (the North Branch, the North Channel, the North Passage and the South Passage) into the sea. In order to meet the needs of large vessels getting in and out of the Yangtze Estuary, the Yangtze Estuary Deepwater Channel (YEDC) Regulation Project was implemented by a combination of training works and dredging in three phases. From 1998 to 2005, the Projects of Phase I and Phase II, including the construction of the diversion project of the South Passage and North Passage, South Jetty and North Jetty with the lengths of 48 km and 49 km respectively and 19 spur dikes with a total length of about 30 km, were completed. The layout of the main training building is shown in Fig. 2. The depth of the North Passage reaches 8.9 m after the completion of Phase I and 10 m after the completion of Phase II. Dredging in Phase III began in September 2006. Due to high backsilting rate and difficulties of deepening the channel, such measurements as the extension of the lengths of spur dikes and the increase of the number of dredging vessels were adopted. The YEDC ( Fig. 1) with the length of 90 km, width of 350-400 m and depth of 12.5 m (below mean sea level) was finally constructed in March 2010. Under the comprehensive effects of the runoff, tide currents, wind waves and salt, the channel backsilting quantity is large and has the following features: (1) the backsilting quantities from June to November account for 80% of the annual backsilting quantity; (2) the backsilting quantity under normal hydrological conditions accounts for 80% of the total backsilting quantity; the backsilting quantities caused by typhoon storm tides account for 20%; (3) the backsilting quantities in the middle and lower sections of the YEDC account for 80% of the whole channel's backsilting quantity. As for the backsilting problems of the YEDC, many scholars studied channel sediment sources, backsilting forecast and siltation reduction measures. Liu et al. (2011) studied the influences of overtopping water and sand of jetties on the channel backsilting and con-sidered that the sediment on the Jiuduansha Shallow, flowing over the South Jetty and entering into the North Passage, made a significant contribution to the backsilting. Based on the data measured near the bottom of the channel, Gu et al. (2017) built a backsilting quantity numerical model of the YEDC and drew a conclusion that the backsilting quantities mainly occurred during the middle and neap tides and less backsilting quantities during the spring tides because of strong hydrodynamic power, significant erosion, and high sediment concentration near the bottom.
Fractal interpolation method, an important content in the fractal theory, can well depict the rough, complex, self similar or self affine objects in real world. It can not only be used to fit smooth curve and stationary data but also show its unique superiorities in the fitting of non-smooth curve and non-stationary data. Therefore, the fractal interpolation method has been widely applied in many fields. Lima et al. (2003) applied IDW interpolation method and a fractal filtering technique to analyze the natural and human factors of the sediment samples in the water system of Campania region (Italy) and the application results are quite acceptable. Liu et al. (2008) used fractal interpolation to predict the monthly runoff of the Yele Reservoir of the Nanya Rriver and the verification analysis shows that the prediction results of the model are good and meet the requirements in practical application. Luor (2017) carried out fractal interpolation analysis for random data sets, showing the ability of fractal method in data analysis. Fractal interpolation method has been widely used in hydrologic data statistics and analyses, but seldom applied in the sediment siltation research.
Because the backsilting of the YEDC is related to dynamic conditions, dredging methods and other factors, the existing theories and simulation methods can hardly describe this phenomenon. This paper, based on the accumulated monthly measured data of the backsilting quantities in the channel, adopts the fractal interpolation method to analyze annual backsilting in the YEDC, providing an approximation method for prediction of the annual backsilting quantity in the YEDC.

Applicability of the fractal interpolation
The precondition that the time series data can be studied by the fractal theory is that the series are fractal. Therefore, it is necessary to analyze the mean, standard deviation, skewness and kurtosis of time series data, as well as the difference relative to the normal distribution.
2.1 Calculation of the characteristic values of the time series The mean of the time series data: x n where, i is the time series, i=1, 2, …, n; is the mean value corresponding the time series data, and x i is the data value for i. The standard deviation of the time series data is: x where, is the mean value.ŝ the skewness of the time series data is the deviation degree of data distribution relative to the normal distribution, it can be calculated as follows: (3) When , the distribution pattern is the same to the normal distribution; the larger indicates the larger deviation from the normal distribution.k The kurtosis of the time series data is the steepness degree of the data distribution, which can be expressed as: In Eq. (4), if , the steepness of the data distribution is the same to that of the normal distribution; If , the data distribution is steeper, or gentler than the normal distribution.

Characteristic values of the backsilting quantity in the YEDC
The main dynamic factors in the Yangtze Estuary are runoff, tide and wave. The runoff is regulated by the Three Gorges Project; although there is a flood or dry season within one year, the discharge varies a little in years. The offshore tides change periodically. The wave is affected by monsoons and typhoons and is periodic. Fig. 1 shows the dredged units of YEDC. According to the monthly dredging volume per dredged unit and the water depth, the monthly backsilting quantities in each unit can be calculated (Fig. 2); and then the monthly backsilting quantities in the whole channel can be obtained. In Fig. 2, the positive value is the erosion quantity and negative value is the siltation quantity. Therefore, the backsilting quantities of two months were equally divided to be used as the backsilting quantities of the relevant months. Eqs.
(1)-(4) are used to calculate characteristic values such as the means, standard deviations, skewness and kurtosis of monthly backsilting quantities in the YEDC. The skewness and kurtosis of the normal distribution are 0 and 3, respectively. From the statistics on the characteristic values (seen in Table 2) of annual backsilting quantities, the mean value of annual backsilting quantity in 2012 is the maximum; and then the means of annual backsilting quantities is decreasing and approaching. The maximum standard deviation is 737.66, and the minimum is 457.51; compared with 2016-2017, the distribution of monthly backsilting quantity from 2011 to 2015 is more dispersed. The skewnesses of 2011 and 2013-2015 are smaller than 0; and compared with the normal distribution, they are left skewed; the skewnesses of 2012 and 2016 are larger than 0, showing that they are right skewed. Except that the kurtosis in 2013 is larger than 0 and its distribution is even steeper than the normal distribution, the kurtosises of other years are smaller than 0 and the distributions are even smoother than the normal distribution. The time series of monthly backsilting quantities in the channel from 2011 to 2017 are not consistent with the normal distribution and have fractal characteristics so that they can be studied by using the fractal theory.

Fractal interpolation method
Traditional mathematical interpolation fitting function is mostly expressed by such elementary functional combinations as polynomial, rational function or trigonometric func-tion. However, the fractal interpolation function is implemented by using iteration function system (IFS). The fractal interpolation can be sufficiently close to 1-2 or 2-3 dimensional images; and its fractal dimension approaches to the fractal dimensions of these data within the appropriate scale.
As for the fractal interpolation, an iterative function system should be constructed so that the target point set A is a part of the interpolation function. The specific method, with a target point set A being constructed and equal to the iterative function system of the interpolating function f(x), is as follows. { The observation data set is known, where, x i <x i +1, the geometric figure of the interpolation function f(x) corresponding to this data set passes through each point continuously, that is, y i =f(x i ). Barnsley (1986) gave an interpolation function f(x) constructed by using IFS; ω i in IFS has the following form of the affine transformation: That is: Eq. (6) is to ensure that the left endpoint (x 0 , y 0 ) of the large interval is mapped to the left endpoint of the sub-interval (x i-1 , y i-1 ) and the right endpoint (x N , y N ) of the large interval is mapped to the right endpoint of the sub-interval (x i , y i ). Parameter d i is the vertical scaling factor of the affine transformation; and its size affects the fluctuation degree of the function values. Let d i be the free parameter (d i <1, otherwise, the IFS does not converge) and Eq. (5) substitute into Eq. (6), and then the parameters can be expressed as: When the interval is determined, a i and e i can be obtained respectively from Eqs. (7) and (8). There are three unknowns in Eqs. (9) and (10), so the iterative calculation is required.
For simplicity, a random factor is adopted to determine   .
where, y max and y min are the interval maximum and minimum values obtained by extending n 0 points forward and backward with the pending interpolation point as the center.
ε=1+rand (1), and rand (1) is a random number between 0 and 1 generated by the computer. When d i is calculated from Eq. (11), c i and f i can be solved by Eqs. (9) and (10). After the above parameters being obtained, x and y on the right side of Eq. (5) are assigned to the initial values respectively. After multiple iterations, ω i (x) and ω i (y) of the interpolation points can be obtained. As the number of iterations increases, the fitting degree of the interpolation curve and original curve increases unceasingly; after many iterations, a stable interpolation curve is formed. The curve not only passes through the interpolation points and also greatly approaches to the original curve (Wang and Li, 2008).
When the fractal interpolation method is used, the measured monthly backsilting quantity is used as a control node. The adjacent monthly backsilting's interpolation curve presents "rough, complex, self-similar" characteristics of the control points. It is not consistent but similar to the time series distribution characteristics of actual monthly backsilting quantity, increasing the fractal dimensions among monthly backsilting quantities.

Influences of the vertical scaling factor
Take the backsilting quantities in the channel from January 2011 to December 2017 as the study objects; calculate the corresponding monthly parameters a i , b i , c i , d i , e i and f i in accordance with Eqs. (7)-(11); take the monthly backsilting quantities in January and December of each year as the boundary values x 0 and y 0 as well as x N and y N and then substitute into Eq. (5). The fractal interpolation is obtained by three iterations of ω i (x) and ω i (y).
The backsilting quantity Q can be calculated by integration of the backsilting quantities obtained by the fractal interpolation during time t: where, Q i is the backsilting quantity obtained by the i-th interpolation, dt i is the time interval of the i-th interpolation and M is the number of the measured points. When the fractal interpolation is performed for each month during seven years, M=12×7=84; when the fractal interpolation is performed for 12 months, M=12. With Eq. (12), the fitting line of backsilting quantities from 2011 to 2017 is obtained. Fig. 3 shows the fitting line of the backsilting quantities in the YEDC in 2014.
(1) Relationship between d i and monthly backsilting quantities Take the monthly backsilting quantities of the YEDC in 2014 to conduct the study. Use d i =-0.5, -0.3, 0, 0.3 and 0.5 to respectively calculate the correspondingly annual backsilting quantities. From Fig. 4, for d i >0, the backsilting quantities increase with the increasing d i ; for d i <0, the backsilting quantities decrease with the increasing d i , all fluctuating around the measured values.
For the analysis of the relationship between monthly backsilting quantities and d i from 2011 to 2017, taking the monthly backsilting quantities in each of the January from 2011 to 2017 as an example, seven monthly backsilting quantities in January are used to conduct the fractal interpolation; after three iterations, 6 3 (=216) backsilting quantities are obtained. Taking the average value after the integration, the mean backsilting quantity in January is obtained. This backsilting quantity is obtained based on the fractal interpolation, which is slightly different from the measured one. The backsilting quantities in other months can be calculated  in the same way. The corresponding relationship between the average monthly backsilting quantity and d i is shown in Fig. 5. It can be seen that there is a linear relationship between the monthly backsilting quantities and d i . The line types of January and April, February and March are similar respectively. The slopes of May, June, November and December are close, whereas the slopes of July, August, September and October are close. The influences of d i on the monthly backsilting quantities are reflected in the intercept and slope.
(2) Relationship between d i and the annual backsilting quantity Analysis on the relationship between the annual backsilting quantity and d i from 2011 to 2017 is carried out with d i =-0.5, -0.4, -0.3, -0.2, -0.1, 0, 0.1, 0.2, 0.3, 0.4 and 0.5. Taking year 2011 as an example, the monthly backsilting quantities from January to December are used to make the fractal interpolation; after three iterations, 11 3 (=1331) data of the backsilting quantities are obtained. The whole backsilting quantity in 2011 is obtained by the integration.The annual backsilting quantity obtained from the fractal interpolation is slightly different from the measured one. The annual backsilting quantity for rest years can be calculated in the same manner. The corresponding relationship between the annual backsilting quantity and d i can be obtained, as shown in Fig. 6. And then, the fractal interpolation is made to obtain the backsilting quantities from January 2011 to December 2017; the data of 83 3 (=571787) backsilting quantities can be obtained after three iterations. Namely, the corresponding relationship between the average annual backsilting quantity and d i can be obtained, and there is a linear relationship between the backsilting quantity and d i . Therefore, d i can be regarded as a comprehensive impact factor reflecting the backsilting quantity in the YEDC, involving hydrodynamic forces (runoff, tide and wave), silt characteristics (size, settling velocity and viscosity), and dredging operations. Different vertical scaling factors lead to different backsilting quantities.

Prediction of the backsilting quantity in the YEDC
It is very difficult to accurately predict the backsilting quantity in the YEDC, not only due to the changes of offshore dynamic conditions, the water and sediment coming from the upstream, but also related to the dredging frequency and operation mode. This paper aims to predict the backsilting quantity in the coming year on the basis of the collected backsilting quantity data.

Prediction of the annual backsilting quantity
From Fig. 6, the linear relationship between the average siltation quantity from 2011 to 2017 and d i can be obtained. Namely, where, Q is the annual backsilting quantity, and d i is the vertical scaling factor. The vertical scaling factors corresponding to the measured backsilting quantities from 2011 to 2017 respectively are 0.148, 0.353, 0.089, 0.024, -0.067, -0.274, and -0.266, ranging from -0.266 to 0.353. According to the time series correlation of the fractal theory, let where, d N+1 is the vertical scaling factor for the predicted year, and N is the number of the years in which backsilting quantity is known.
Considering the changing trend of the backsilting quantities of the YEDC in recent years, let the mean of vertical scaling factors from 2015 to 2017 to be the vertical scaling factor for 2018, i .e. d i2018 =-0.202; and then, the backsilting quantity in 2018 can be calculated as -5.946×10 8 m 3 by Eq. (14).

Prediction of monthly backsilting quantities
Based on the predicted backsilting quantity in 2018, in  LI Lan-xi et al. China Ocean Eng., 2018, Vol. 32, No. 3, P. 341-346 accordance with the relation between the backsilting quantity and d i shown in Fig. 5, the corresponding monthly backsilting quantities in 2018 (seen in Fig. 7) can be obtained. The backsilting quantities from January to April account for about 2% of the annual backsilting quantity, 15% for May to June, 65% for July to October 65%, and 18% for November to December. It can be seen that the backsilting quantities from June to November account for 80% of the annual backsilting quantity.

Conclusions
This paper adopts the fractal interpolation method to study the measured backsilting quantities in the YEDC and the following conclusions can be drawn.
(1) According to the data of the monthly backsilting quantity for 7 years measured at the YEDC, by analyzing the characteristic values such as the means, standard deviations, skewness and kurtosis, it can be concluded that the fractal interpolation method is applicable to study the temporal distribution of backsilting quantity.
(2) Through the fractal interpolation on the backsilting quantities from January 2011 to December 2017, when the vertical scaling factor d i is larger than 0, the backsilting quantity increases with d i ; when the vertical scaling factor d i is smaller than 0, the backsilting quantity decreases with d i .
(3) Through the fractal interpolation, the annual siltation quantities and monthly siltation quantities respectively have a linear relationship with the vertical scaling factor d i and are directly proportional to the slope and intercept.
(4) Based on the fractal interpolation of the measured backsilting quantities from 2011 to 2017, a calculation formula for prediction of the backsilting quantity in the YEDC is obtained, which can be used to predict the annual backsilting quantity and monthly backsilting quantities of the channel in future years.
(5) Since the YEDC has been run only for 7 years, the calculation formula for prediction of the backsilting quantity still needs to be verified by the backsilting quantities in future.