Effect of Trends on Detrended Fluctuation Analysis of Precipitation Series

Abstract

We use detrended fluctuation analysis (DFA) method to detect the long-range correlation and scaling properties of daily precipitation series of Beijing from 1973 to 2004 before and after adding diverse trends to the original series. The correlation and scaling properties of the original series are difficult to analyze due to existing crossovers. The effects of the coefficient and the power of the added trends on the scaling exponents and crossovers of the series are tested. A crossover is found to be independent of the added trends, which arises from the intrinsic periodic trend of the precipitation series. However, another crossover caused by the multifractal vanishes with the increasing power of added trends.

1. Introduction

Many physical and biological systems exhibit complex behavior characterized by long-range power-law correlations. Traditional approaches such as the power-spectrum and scaled-Hurst analysis are limited to quantify correlations in stationary signals. In recent years, detrended fluctuation analysis (DFA) has been established as an important tool for the detection of long-range correlations in time series with nonstationarities. DFA is a scaling analysis method providing a quantitative parameter, the scaling exponent 𝛼, to represent the long-range autocorrelation properties of a signal. The advantages of DFA over many other methods are that it permits the detection of correlations in apparent nonstationary time series and also avoids the spurious detection of seemingly long-range correlations that are artifact of nonstationarity. DFA which is a nonparametric approach for data mining has been successfully applied to diverse fields of interest such as DNA, heart rate dynamics, neuron spiking, human gait, cloud structure, economical time series, and long-time weather records as well as [1–9]. Besides, many parameter models as well as relevant prediction also have been systematically explored such as in traffic flows with remarkable results [10–14].

A fact exists that precipitation has a dramatic effect on agriculture and plays a significant role in human’s activities. The study of precipitation can be utilized for several purposes, including hydrological structure design, flood prevention, and so forth. Precipitation has been long analyzed by traditional statistics, and effective methods as well as prediction models have been developed in bulk to investigate its role [15–17]. There also exist many investigations of scaling behaviors and multifractal characterization of the precipitation records [18–20]. However, traditional time series analysis of precipitation always produces spurious results due to the highly nonstationary nature of precipitation signals. Matsoukas et al. [21] used detrended fluctuation analysis to quantify the correlation properties of precipitation time series but did not describe them in detail.

In the paper, we detect the long-range correlations of the daily precipitation series collected from 21 weather stations of Beijing through about 30 years and investigate their correlation properties together with the influence of added trends under the method of DFA. As external trends are the main components which affect the correlation properties of a time series, people are trying to eliminate them to gain proper insight into the records. However, in most cases it is difficult to distinguish the trends from the intrinsic fluctuations in data. We add diverse trends on the contrary to the original data and systematically analyze their effect on the correlation properties. The essence of adding the trends in the paper is a preprocessing as the trends will be the functions of original series.

The organization of this paper is as follows: in Section 2, we briefly introduce the DFA method. Section 3 is about the details of the precipitation data we used in this paper. In Section 4 we detect the correlation properties by calculating the scaling exponents and crossover times of the original series before and after adding correlated trends. We summarize in Section 5.

2. Methodology

Experimental series are often affected by nonstationarity and fractality [22]. To investigate the scaling behavior of fluctuations, external trends are expected to be well distinguished from the intrinsic fluctuations of the system. If trends exist in the data, Hurst rescaled-range analysis and other nondetrending methods might give spurious results [5–8]. Very often we never know the reasons for underlying trends in collected data and even worse the scales of the underlying trends. DFA is a well-established and robust method for determining the scaling behavior of noisy data in the presence of diverse trends [17–21].

For a record {𝑋(𝑖)},𝑖=1,2,…,𝑁, where 𝑁 denotes the length of record, the DFA procedure briefly involves the following four steps.

Step 1 1. We determine the profile 𝑥(𝑖), 𝑥(𝑖)=^𝑖_𝑘=1𝑋_𝑘−⟨𝑋⟩, (2.1) where ⟨𝑋⟩ is the mean of the record.

Step 2. We cut the profile 𝑥(𝑖) into 𝑁_𝑠=[𝑁/𝑠] boxes of the same size 𝑠. In each box, we fit the integrated time series by using a polynomial function, 𝑝_𝑣(𝑖), which is regarded as the local trend. For order-𝑛 DFA, 𝑛 order polynomial function is applied of the fitting approximation. We subtract the local trend in each box and get the detrended fluctuation function 𝑥_𝑠(𝑖): 𝑥_𝑠(𝑖)=𝑥(𝑖)−𝑝_𝑣(𝑖). (2.2)

Step 3. In each box of size 𝑠, we calculate the root mean square (rms) fluctuation 𝐹(𝑠): 𝐹(𝑠)=1𝑁_𝑠^𝑁𝑠_𝑖=1𝑥²_𝑠(𝑖). (2.3)

Step 4. We repeat this procedure for different box sizes 𝑠 (different scales).
If a power-law relation exists between 𝐹(𝑠) and 𝑠, 𝐹(𝑠)∼𝑠^𝛼. (2.4) It indicates the presence of scaling property. The parameter 𝛼, called the scaling exponent or fluctuation exponent, represents the correlation properties of the data. For correlation exponent 𝛾, which is derived from the autocorrelation function, a similar approximation for 𝐹(𝑠) is 𝐹(𝑠)∼𝑠^1−𝛾/2. (2.5) Comparing with (2.4), we find 𝛾=2−2𝛼 for 0<𝛾<1. A brief certification of the relation of 𝛾 and 𝛼 is proposed in [23]. We can determine the correlation exponent 𝛾 by measuring the fluctuation exponent 𝛼. If 𝛼=0.5, there is no correlation (white noise); if 𝛼<0.5, the data is anticorrelated; if 𝛼>0.5, the data is long-range correlated.

3. Data Description

The precipitation data here is collected from 21 weather stations of Beijing from January1 1973 to December 31 2004, 11688 days, as illustrated in Figure 1. A vision processing based on coded structure light has been investigated to acquire 3D data which can be referred for further analysis if necessary [24, 25]. There may not be any precipitation record which is different from such as temperature series. There may be a little precipitation that it is not necessary for the weather stations to record in detail but adopting a word “minim’’ instead of a specific quantity. For convenience, we regard the quantity of the days without precipitation as 0 and the “minim’’ as 0.5. We treat the mean value of the 21 records every day of different stations as a new precipitation series for analysis.

Figure 1 

illustrates the precipitation data of every day in the record years. The unit of the precipitation is 0.1 mm. From this figure, we can clearly see the seasonal trend of the data. We derive similar results from the DFA of the original data shown in next section.

4. Data Analysis

4.1. DFA of the Original Series

First, we detect the correlation behavior of the original series. To get more information, we use the DFA arranging from 1st to 5th order. The original series is a multifractal according to Figure 2(a) since the scaling exponents of each order-𝑛 DFA change twice, that is, two crossovers. At small scales 𝑠, the deviations grow stronger with the increasing DFA order 𝑛. To decrease the impact of the deviations on the calculating of scaling exponents 𝛼^(𝑛), we ignore some small 𝑠 while fitting the curve.

(a)

(b)

Figure 2 

(a) DFA of the original series. For each order-$n$ DFA there exists two crossovers $s_{1x}^{(n)}$ and $s_{2x}^{(n)}$ (We use$s_x^{(n)}$ for them) which divide the curves into 3 different scaling segments whose scaling exponents are: ${\alpha }_0^{(n)}\text{,}$${\alpha }_1^{(n)}$ and ${\alpha }_2^{(n)}$ respectively (We use ${\alpha }^{(n)}$ for them). (b) The DFA of a sinusoidal series given by the function $10*\hspace{0.17em}\sin (20\pi i/N)$, where $N$ is the length of the original series

Crossover times 𝑠^(𝑛)_1𝑥 and 𝑠^(𝑛)_2𝑥 are determined by the intersection of linear fits done on both sides of the crossovers. We choose the point at the intersection in scales {𝑠} as the crossover times 𝑠_𝑥^(𝑛), calculate the slope on both sides of 𝑠_𝑥^(𝑛) to get 𝛼^(𝑛), and exhibit 𝑠_𝑥^(𝑛), 𝛼^(𝑛) in Table 1.

Table 1 

The crossover times 𝑠_𝑥^(𝑛) and the scaling exponent 𝛼^(𝑛) of the original data.

𝑠^(𝑛)_1𝑥 and 𝑠^(𝑛)_2𝑥 of each order-𝑛 DFA divide the series into three different scaling segments. At the first segment where the time scales 𝑠<𝑠^(𝑛)_1𝑥, 0.5<𝛼₀^(𝑛)<1 corresponds to a long-range correlation behavior which indicates that a relatively large magnitude is likely to be followed by a large magnitude event. For the two segments on the two sides of the second crossover 𝑠^(𝑛)_2𝑥, scaling exponents 𝛼^(𝑛) changes from 𝛼₀^(𝑛)>1 to 𝛼₂^(𝑛)<0.5 which means that at large time scales 𝑠, the series is anticorrelated. The scaling behavior of the original series at large scales 𝑠 is similar to the DFA of a sinusoidal series shown in Figure 2(b). Comparing Figure 2(a) with Figure 2(b), there exist some common properties at large scales; after a significant crossover, both scaling exponents turn rather small. It manifests that periodic trend dominates the scaling property at large scales after 𝑠^(𝑛)_1𝑥 which is accordant to the investigations in [4]. It also can be referred that 𝑠^(𝑛)_2𝑥 is highly possible to be dominated by the seasonal trend in the precipitation series.

4.2. DFA of the Series with Correlated Trend 𝐴[𝑋(𝑖)]^𝑝

The complex properties of the original series result from the crossover times 𝑠^(𝑛)_1𝑥 and 𝑠^(𝑛)_2𝑥, so it is difficult to understand the scaling behavior and make a valid prediction. A crossover usually can arise from a change in the correlation properties of the series at different time scales, that is, multifractal, or can often arise from external trends in the data [26]. Diverse methods provide inspiration to produce discrete sequences and continuous functions [27, 28] for simulation. In most cases, people generate long-range correlated experimental data with modified Fourier filtering [29] or “ARFIMA’’ [30] method and superimpose diverse trends on them which are the function of time, like linear, sinusoidal, and power-law trend. The trends effects on the original series are tested, mainly including the crossover, complete with diverse detrending methods based on such as SVD [9], EMD [31–33], Fourier-DFA [34], wavelet analysis [35], and “superposition rule” [4]. In real data, the type of trend is analyzed and proper detrending method is employed correspondingly. It is an attractive and logical direction in solving the crossover caused by trends and deriving a constant scaling exponent. However, at times the type of the trend is difficult to identify, and we note that the information of the series is not fully uncovered just by the DFA method. Here we firstly propose a new method to preprocess the data by adding a trend which is a function of the original series. Then we test whether and how correlated trends added to the original series will affect the correlation properties. A common power-law function 𝐴[𝑋(𝑖)]^𝑝 is used, where 𝐴 is a coefficient and 𝑝 presents the power. We apply DFA method to the new series 𝑌(𝑖)=𝑋(𝑖)+𝐴[𝑋(𝑖)]^𝑝. It is apparent that 𝑌 share the same period of trend with 𝑋 but with different fluctuation magnitudes. As we will see later, 𝑠^(𝑛)_2𝑥 dominated by periodic trend is independent of added trend 𝐴[𝑋(𝑖)]^𝑝, but 𝑠^(𝑛)_1𝑥 rises from the multifractal will disappear. We make 𝐴 and 𝑝 variables, respectively, to operate our study as follow.

4.2.1. Effect of Power 𝑝 on DFA of 𝑌(𝑖)

In this section, 𝑝 is a variable and 𝐴 is a constant 1.

(1) 𝑝 Is Positive Integer
Considering the capability of order-𝑛 DFA in removing trend of (𝑛−1)th order, we give 𝑝 integer values ranging from 1 to 6 to take a whole view of effects of DFA on the series 𝑌(𝑖).
For 𝑠<10³ every order-𝑛 DFA in Figures 3(a) and 3(b), two crossovers exist while in (c) the crossovers only remain in DFA3, DFA4, and DFA5 together with their positions being much closer. For (d), (e), and (f), although the crossovers still exist, one cannot identify them without checking carefully. For each order-𝑛 DFA in (d), (e), and (f), only one crossover exists, which is still marked by only 𝑠^(𝑛)_2𝑥 for convenience. We illustrate the representative crossovers by arrows. It seems that crossovers 𝑠^(𝑛)_2𝑥 in all six subfigures share an identical position. We apply the method in Section 4.1, calculate these crossovers and scaling exponents, and specify the crossover time scales in Table 2. Since the scaling behavior of DFA3, DFA4, and DFA5 in Figure 3(c) on both sides of 𝑠^(𝑛)_1𝑥 is just the same, we calculate one scaling exponent 𝛼₁^(𝑛) for each of them before 𝑠^(𝑛)_2𝑥 in Table 2.
The positions of crossover times 𝑠^(𝑛)_2𝑥 after adding diverse trends are identical to that of the original series, which demonstrate that 𝑠^(𝑛)_2𝑥 are independent of the power 𝑝 of the correlated trend. With the increasing of 𝑝 from 1 to 6, the values of all scaling exponents 𝛼₀^(𝑛) and 𝛼₁^(𝑛) before 𝑠^(𝑛)_2𝑥 decrease to 0.5 while all values of 𝛼₂^(𝑛) after 𝑠^(𝑛)_2𝑥 have the trend to increase to 0.5. To get a clearer view of this phenomenon, for 𝑝=1 and 𝑝=2, we calculate a new scaling exponent 𝛼₁^(𝑛) before 𝑠^(𝑛)_2𝑥 of each order-𝑛 DFA by linear fit. They are specified, respectively in Table 3.
The decreasing trend of 𝛼₁^(𝑛) to 0.5 and the increasing trend of 𝛼₂^(𝑛) to 0.5 are shown in Table 3. With the increasing of 𝑝 the approaching pace of 𝛼₂^(𝑛) to 0.5 seems to be slower than that of 𝛼₁^(𝑛).
𝑠^(𝑛)_1𝑥 vanishes with increasing 𝑝 of the trend 𝐴[𝑋(𝑖)]^𝑝 but 𝑠^(𝑛)_2𝑥 is independent of the added trend. It agrees well with the previous analysis that 𝑠^(𝑛)_1𝑥 and 𝑠^(𝑛)_2𝑥, respectively arise from the multifractal and periodic trend in precipitation series.
The scaling behavior of 𝑌(𝑖) is similar to that of original signal when 𝑝=1, as 𝑌(𝑖)=(𝐴+1)𝑋(𝑖). When 𝑝>1, since the magnitude of 𝑋(𝑖) is mostly larger than 1, the scaling properties of 𝐴[𝑋(𝑖)]^𝑝 will play a vital role and the influence of 𝑋(𝑖) can be negligible which also can be inferred from the scale of the fluctuation 𝐹_𝑛(𝑠). In fact, Figures 3(b)–3(f) demonstrate the fluctuations of 𝐴[𝑋(𝑖)]^𝑝 as well. We note that the scaling properties of 𝐴[𝑋(𝑖)]^𝑝 is also attractive, as it is the function of original series which can be treated a preprocessing.

Table 2 

Scaling exponents 𝛼^(𝑛) and crossover times 𝑠_𝑥^(𝑛) of the series with correlated trends 𝐴[𝑋(𝑖)]^𝑝, 𝐴=1 and 𝑝 are integers from 1 to 6.

Table 3 

Comparison of 𝛼₁^(𝑛) and only 𝛼₂^(𝑛) from 𝑝=1 to 𝑝=6 after recalculating of 𝑝=1 and 𝑝=2.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 3 

DFA of the series with correlated trend $A{[X(i)]}^p$, that is, $Y(i)\mathrm{\hspace{0.17em}\hspace{0.17em}}=X(i)+\mathrm{\hspace{0.17em}\hspace{0.17em}}A{[X(i)]}^p$, where $A=1$, (a) $p=1$, (b) $p=2$, (c) $p=3$, (d) $p=4$, (e) $p=5$, (f) $p=6$.

(2) 𝑝 Is Number with Decimal and 𝑝>1
The crossovers 𝑠^(𝑛)_2𝑥 in Figure 3 are independent of power 𝑝 but 𝑠^(𝑛)_1𝑥 are not which is, respectively, affected by the multifractal and periodic trend of original precipitation series. Figures 3(b), 3(c), and 3(d) show a “vanishing’’ process of 𝑠^(𝑛)_1𝑥, where 𝑠^(𝑛)_1𝑥 remains by adding an order-𝑝 trend while disappears by an order-𝑝+1 trend. However, we also expect to track the whole process of vanishing with a decimal of 𝑝. We test it in Figure 4. First we find some common properties in three pair figures that these significant processes all take place between 𝑚+0.8 to 𝑚+0.9, 𝑚=1,2,3, which are very close to 𝑚+1. As we use DFA to 5th order, the “vanishing’’ process of 𝑠^(𝑛)_1𝑥 just completes when the power of the trend 𝐴[𝑋(𝑖)]^𝑝, 𝑝, is 3.9, close to 5−1=4. Thus we can make a hypothesis that if our DFA is order 𝑛+1, there will be some significant process around each 𝑝=𝑖+0.9, 𝑖=1,2,…,𝑛−1, close to 𝑖+1 and we get the final stable property around 𝑝=𝑛−0.1, close to 𝑛. It may be related to the fact that an (𝑛+1)th DFA can eliminate an 𝑛th order polynomial trend due to the integration in DFA algorithm.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 4 

DFA of the series with correlated trend $A{[X(i)]}^p,$ that is, $Y(i)=X(i)+\mathrm{\hspace{0.17em}\hspace{0.17em}}A{[X(i)]}^p,$ where $A=1$, (a) $p=1.8$. (b) $p=1.9$. (c) $p=2.8$. (d) $p=2.9$. (e) $p=3.8$. (f) $p=3.9$. (a)$\to$(b) is the process that crossover $s_{1x}^{(1)}$ changes from $84(<100)$ to $106(>100)$. (c)$\to$(d) is the “vanishing’’ process of $s_{1x}^{(1)}$ and $s_{1x}^{(2)}$ while (e)$\to$(f) is the “vanishing’’ process of $s_{1x}^{(n)},$$n=\mathrm{3,4},5$. Some typical crossovers are tagged by arrows.

(3) 𝑝<1
As we can see from Figure 3(a), when 𝑝=1, that the extra trend 𝐴[𝑋(𝑖)]^𝑝 does not have any influence on the correlation properties of the series. And as most data in our original precipitation, series are larger than 1.0 (0.1 mm). We can guess that when 𝑝<1, that is, the adding correlated trend is weaker than that of one above, the correlated trend does not affect the correlation properties either and we prove it in Table 4.
We find from Table 4 that the crossovers 𝑠_𝑥^(𝑛) of 𝑝<1 are just the same as the original data. They are rather close to the original data as well by observing their scaling exponents. The trends are so weak that there is little influence of 𝐴[𝑋(𝑖)]^𝑝 on the correlation properties of 𝑋(𝑖).

Table 4 

DFA of the series with correlated trend 𝐴[𝑋(𝑖)]^𝑝, that is, 𝑌(𝑖)=𝑋(𝑖)+𝐴[𝑋(𝑖)]^𝑝𝑝<1. We pick four values of 𝑝 as below.

4.2.2. Effect of coefficient A on DFA of 𝑌(𝑖)

In Section 4.2.1, we have studied the correlation properties of the series with correlated trend 𝐴[𝑋(𝑖)]^𝑝, where 𝐴 is a constant 1. And we find that the crossovers 𝑠^(𝑛)_2𝑥 are independent of power 𝑝 of trend 𝐴[𝑋(𝑖)]^𝑝. So how about coefficient 𝐴?

In this section we vary 𝐴 to find the relation between 𝐴 and the correlation properties of 𝑌(𝑖). As the correlation properties for DFA of 𝑌(𝑖) in Section 4.2.1 is the most complicated when 𝑝=2, here we just test coefficient 𝐴 with 𝑝=2.

There is a clear vision that the crossover times of (a), (b), (c), and (d) are just the same as Figure 3(b). And the scaling exponents 𝛼^(𝑛) seem to be also the same that we exhibit them in Table 5 to make it evident. The scaling exponents 𝛼^(𝑛) of Figures 5(a), 5(b), 5(c), and 5(d) are not exactly the same by values but they hold the same scaling behavior with Figure 3(b).When 𝐴 is large, the scaling exponents 𝛼^(𝑛) are rather close to the ones of 𝐴=1 in Figure 3(b). The crossovers in DFA result from the competition between the scaling of original series and the scaling of the trend. The case of 𝐴<0 is similar to that of 𝐴>0 which can be inferred by the following analysis. If the original series dominates the scaling behavior, apparently the symbol of 𝐴 can be negligible. When the added trends prevail, most values may change the signs, but the magnitude of fluctuation and the period of trend remain unchangeable. Of course, 𝐴=−1 and 𝑝=1 should be excluded. So the correlation properties with correlated trend 𝐴[𝑋(𝑖)]^𝑝 are independent of the coefficient 𝐴.

Table 5 

Scaling exponents of Figures 5(a), 5(b), 5(c), and 5(d).

(a)

(b)

(c)

(d)

Figure 5 

DFA of the series with correlated trend $A{[X(i)]}^p,$ that is, $Y(i)\mathrm{\hspace{0.17em}\hspace{0.17em}}=X(i)+\mathrm{\hspace{0.17em}\hspace{0.17em}}A{[X(i)]}^p$ where $p=2$, (a) $A=0.01$, (b) $A=10$, (c) $A=100$, (d) $A=10000$.

4.3. DFA of the Series with Other Correlated Trends

Discussion of Section 2 has told that crossovers 𝑠^(𝑛)_2𝑥 of the original data are independent of the trend 𝐴[𝑋(𝑖)]^𝑝 which presents another question which is are they still independent of other correlated trends? In this section we test another common trend 𝐴𝑒^{𝑝𝑋(𝑖)/𝑥} (𝑥 is the max of 𝑋(𝑖)) which is stronger than 𝐴[𝑋(𝑖)]^𝑝 with the same method applied in Section 2. We test the effects of 𝐴 and 𝑝 on the correlation properties of 𝑍(𝑖)=𝑋(𝑖)+𝐴𝑒^{𝑝𝑋(𝑖)/𝑥}, respectively, in Figures 6 and 7.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 6 

DFA of the series with correlated trend $A{e}^{pX(i)/\overline{x}}$, where $A=10$, (a) $p=2,$ (b) $p=3,$ (c) $p=5,$ (d) $p=6,$ (e) $p=7,$ (f) $p=8$.

(a)

(b)

(c)

(d)

Figure 7 

DFA of the series with correlated trend $A{e}^{pX(i)/\overline{x}},$ where $p=5$, (a) $A=0.01$, (b) $A=1$, (c) $A=100$, (d) $A=10000$.

Although Figure 7 indicates that the correlation properties of series with correlated trend 𝐴𝑒^{𝑝𝑋(𝑖)/𝑥} are dependent of 𝐴 which are not the same as they are for trend 𝐴[𝑋(𝑖)]^𝑝, we can see form Figures 6 and 7 that the crossover times 𝑠^(𝑛)_2𝑥 are just identical to the ones of original series for each figure. So crossovers 𝑠^(𝑛)_2𝑥 of the original data are independent of the correlated trend 𝐴𝑒^{𝑝𝑋(𝑖)/𝑥}. It is an obvious conclusion that crossovers 𝑠^(𝑛)_2𝑥 are independent of trend 𝐴log(𝑝𝑋(𝑖)) as well, because it is weaker than 𝐴[𝑋(𝑖)]^𝑝, not to mention 𝐴𝑒^{𝑝𝑋(𝑖)/𝑥}. Thus we can see that crossovers 𝑠^(𝑛)_2𝑥 of the original data are rather significant and they are so stable that it is necessary to take further investigations to gain insight into them.

5. Conclusion

In summary, DFA of the original series indicates the complex correlation properties by two obvious crossovers and three different scaling segments. The second crossover 𝑠^(𝑛)_2𝑥 is proved to be independent of the adding correlated types of trend:𝐴log(𝑝𝑋(𝑖)), 𝐴[𝑋(𝑖)]^𝑝, and 𝐴𝑒^{𝑝𝑋(𝑖)/𝑥},while the first crossover 𝑠^(𝑛)_1𝑥 disappears with added trend. They are induced by different reasons while behave similarly if we just analyze original precipitation series. The paper also provides a method to distinguish the multifractal and trend effects on the scaling behavior. With the development of study on correlation, scaling exponents and crossovers will take more significant roles in providing foundation theories for precipitation series predictions based on correlations theories.

Acknowledgments

The financial supports from the funds of The National High Technology Research Development Program of China (863 Program) (2007AA11Z212), the China National Science (60772036) and MEDF (20070004002) are gratefully acknowledged.