Mathematical Problems in Engineering
Volume 2015 (2015), Article ID 174289, 24 pages
http://dx.doi.org/10.1155/2015/174289

Research Article

Maximum Likelihood Estimation of Clock Skew in IEEE 1588 with Fractional Gaussian Noise

Chagai Levy and Monika Pinchas

Department of Electrical and Electronic Engineering, Ariel University, Ariel 40700, Israel

Received 19 October 2014; Revised 6 February 2015; Accepted 23 February 2015

Academic Editor: Frédéric Hamelin

Copyright © 2015 Chagai Levy and Monika Pinchas. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

To support system-wide synchronization accuracy and precision in the sub-microsecond range without using GPS technique, the precise time protocol (PTP) standard IEEE-1588 v2 is chosen. Recently, a new clock skew estimation technique was proposed for the slave based on a dual slave clock method that assumes that the packet delay variation (PDV) in the Ethernet network is a constant delay. However, papers dealing with the Ethernet network have shown that this PDV is a long range dependency (LRD) process which may be modeled as a fractional Gaussian noise (fGn) with Hurst exponent () in the range of . In this paper, we propose a new clock skew estimator based on the maximum likelihood (ML) technique and derive an approximated expression for the Cramer-Rao lower bound (CRLB) both valid for the case where the PDV is modeled as fGn (). Simulation results indicate that our new clock skew method outperforms the dual slave clock approach and that the simulated mean square error (MSE) obtained by our new proposed clock skew estimator approaches asymptotically the developed CRLB.

1. Introduction

Recently, the topic of synchronization in network nodes has attracted significant attention for a number of applications [1–10]. In many cases the global positioning system (GPS) is used to synchronize frequency and time accuracy in the sub-microsecond range. On the other hand, although GPS can offer reliable clock synchronization accuracy, network operators are seeking to minimize the use of GPS within their networks (especially for the indoor places) due to a weak indoor GPS signal [11]. Technologies like synchronous Ethernet (SyncE) [12], network time protocol (NTP) [13], and IEEE 1588v2 (PTP) [14] have been employed to deliver frequency synchronization (GPS, SyncE, NTP, and PTP) and time synchronization (GPS, NTP, and PTP) in Ethernet networks. This provides significant cost savings in network equipment as well as in ongoing installation and maintenance [15].

The SyncE is defined by the ITU [12] as a means of using Ethernet to transfer frequency via the Ethernet physical layer. It recovers frequency from one or more Ethernet signals and delivers timing to outgoing Ethernet signals. However, in order to use SyncE effectively, hardware changes in the switches and routers across the network are needed, which are not yet done.

The NTP was proposed for the first time in 1985 and has been revised several times [16]. The time filtering rules and the servoregulation algorithms implemented by the NTP clients have been designed for large networks; hence they are able to trace time down to a few milliseconds on geographically large-scale network [16]. However, the NTP protocol is not designed for local area networks (LANs) and, due to slow response and software clock implementations (timestamping at application layer), cannot compete with PTP performance even when special implementation are used [17].

IEEE 1588 v2 is also named as PTP, used for time synchronization in networked measurement and control systems. IEEE 1588 is a master-slave time synchronization protocol that targets clock synchronization accuracy of less than 1 μs [18] due to the precise timestamping at the physical layer. In the last decade, the PTP has experienced a large diffusion in systems that require high performance synchronization over Ethernet, since its simplicity and autoconfiguration are suitable for embedded and low cost applications [19].

There are two main functions established by the standard [20]: (i) establishing a master-slave hierarchy of clocks in which each slave synchronizes to its master and (ii) making the necessary information available for the slaves to perform this synchronization. Most of the basic approaches [20, 21, 23] involve clock synchronization that is based on the knowledge of the link propagation delay that is the master-to-slave message delay (uplink) and slave to master message delay (downlink). By assuming symmetric communication links (uplink delay = downlink delay) [20, 21, 23] estimate the clock offset, that is, the time difference between the master and the slave. However, in a realistic network the up- and downlink delays may be very different. Thus, this approach may lead to errors in clock synchronization in a realistic network. Recently, a clock synchronization method (frequency synchronization) was proposed [18] based on dual slave clocks in a slave, suitable for symmetric and asymmetric communication links. This estimation method [18] is based on the knowledge of the one-way delay (OWD) (on the uplink delay) instead on the knowledge of the two-way delay (2WD) (namely, on both directions up- and downlink delays) as introduced in the basic concept of the IEEE 1588 [22]. It does not assume having symmetric communication links (uplink delay = downlink delay) unlike [20, 21, 23]. Thus, its estimation method accuracy [18] does not depend on the up- and downlink symmetry unlike in [20, 21, 23]. Therefore, if the up- and downlink are asymmetric, the estimation method accuracy of [18] will not be affected by this asymmetry unlike how it is in the case of [20, 21, 23]. However, [18] does not take into account the PDV existing in the network. According to [24], PDV is a random delay due to the behavior of switches or routers in the network. The primary source of this random delay is the output queuing delay, caused when a PTP message arrives at a switch or router when the exit port is blocked by other traffic, and the PTP Sync message has to wait in a queue. This random delay (PDV) is characterized according to [25–30] as LRD characteristics such as fGn with Hurst exponent in the range of .

According to [31], fGn is a commonly used model of network traffic with LRD. In [32], the authors state that fGN with a single parameter (with ) is a widely used traditional traffic model and refer the reader to [33–35]. They [32] also add that traditional methods to synthesize traffic are based on fGN with a single parameter and refer the reader to [36–40]. Please note that, in [30, 41], fGN was used to model the network traffic where the work in [26] is a supplement of [41] that supplies experimental investigations to show that fGN can yet be used for modeling autocorrelation functions of various types of network traffic consistently in the sense that the modeling accuracy (expressed by mean square error) is in the order of magnitude of 10⁻³. The conclusion of [30] was that the analyses of generating fGN data have shown that this generator should be recommended for practical simulation studies of high-speed telecommunication networks (Ethernet, ATM, VBR video traffic, Web traffic, Telnet, and FTP), since it is very accurate.

It should be pointed out that there are also other works [42–48] that model the PDV as a Gaussian or exponential delay. These papers [42–48] do not need the symmetric assumption on the up- and down link unlike [20, 21, 23]. However, they ([42–48]) work with the 2WD messaging technique where both sides (uplink and downlink) are affected by the PDV unlike with the OWD messaging technique. Therefore, in cases where the uplink PDV is lower than the downlink PDV we may get with the OWD messaging method which is affected only by the uplink PDV, better performance compared with the 2WD messaging technique. As already explained earlier in this paper, [18] uses the OWD messaging technique (unlike [42–48]) and does not assume a symmetric communication link unlike [20, 21, 23]. But, [18] does not take into account the PDV. Please refer to Table 1 for a brief comparison between the mentioned methods. Our motivation in this paper is to develop a new normalized frequency estimator for the clock skew that works with the OWD messaging technique, does not assume a symmetric communication link, and takes into account the PDV modeled with the fGN.

Table 1: Comparison between the various estimation methods.

Model
Reference
Up- and Downlink assumption
PDV model existence in the algorithm
PDV model
Exposure to PDV (2WD/OWD)

Basic concept of IEEE 1588
[20, 21, 23]
Symmetry
—
—
Up- and downlink

Dual slave clock in a slave
[18]
Asymmetry/symmetry
—
—
Uplink only

Gaussian delay model
[42, 46]
Asymmetry/symmetry
✓
Gaussian
Up- and downlink

Exponential delay model
[42–45, 47, 48]
Asymmetry/symmetry
✓
Exponential
Up- and downlink

MLE of clock skew in IEEE 1588 with fGn
(New model)
Asymmetry/symmetry
✓
fGn
Uplink only

In this paper we propose a new normalized frequency estimator for the clock skew based on the dual slave clocks technique proposed in [18], where the PDV is modeled as fGn unlike how it was done in [18]. The clock skew is defined as the normalized difference in the clock frequencies between the slave and the master. This new method is based on the ML estimation technique and outperforms the existing method [18] from the MSE point of view in a realistic communication link. In addition we derive in this paper an approximated expression for the Cramer-Rao lower bound (CRLB) depending on the Hurst exponent in the range of . Simulation results indicate that the MSE obtained from our new proposed normalized frequency estimator declines and approaches asymptotically our new derived CRLB.

The paper is organized as follows. The system description is first introduced in Section 2. The proposed ML estimator and the CRLB are presented in Section 3. In Section 4 simulation results are presented and the conclusion is given in Section 5.

2. System Description

2.1. IEEE 1588: Basic Concept

In this subsection we describe the basic concept of IEEE 1588 based on [18]. Figure 1 (recalled from [18] and slightly modified) describes the PTP messaging. The master sends Sync messages at a default rate of once every 2 s. The master clock measures the time when the first Sync message is sent. When a slave clock receives the Sync message, it measures and stores the time of reception . The master then sends a Follow_Up message that contains the actual value of the timestamp . At this point the slave calculates the master-to-slave delay [18]:

On the other hand, the master-to-slave delay is the sum of the time difference between the master and the slave clocks (i.e., clock offset) and the master-to-slave message delay [22]:

Figure 1: PTP messaging timing diagram.

In order to estimate the clock offset the slave periodically sends the Delay_Req message and stores the transmission time with a timestamp . When the master receives the Delay_Req message, it sends a Delay_Resp message, which conveys the timestamp of the reception time of the request message. The slave calculates the slave-to-master delay from these timestamps as [18]

On the other hand, the slave-to-master delay is the sum of the negative time difference between the slave and the master clocks (i.e., clock offset) and the slave-to-master message delay [22]:

Assuming that the master-to-slave message delay is equal to the slave-to-master message delay we may write [22]

where is the estimated message delay of the uplink or downlink. By using (1) to (5), the slave estimates and the clock offset () as

The slave adjusts its time to minimize the clock offset (), thereby synchronizing with the master clock. Please note that (7) assumes that there is no frequency offset between the master and the slave clocks.

Every oscillator has its unique clock frequency consisting of a nominal frequency in addition to a normalized frequency difference of parts per million (PPM) from the nominal frequency multiplied by this nominal frequency. In addition, every oscillator has a random jitter . Therefore the master and the slave clocks may have different frequencies and thus are not synchronized. For example, suppose that mater's frequency is () and slave's frequency is (); thus their frequencies are not synchronized. In that case, the clock offset () generally keeps increasing and may cause the communication link to fail. Thus, the importance of estimating the normalized frequency difference between the slave and the master that is clock skew ().

2.2. Dual Slave Clocks Method

In this subsection we present briefly the dual slave clock method in a slave as presented in [18]. The basic PTP concept makes the assumption that the up- and downlink delays are equal. However, in a realistic network the up- and downlink delays may be very different. Thus, this approach (where the up- and downlink delays are assumed to be equal) may lead to errors in the clock synchronization task in a realistic network. To eliminate the assumption of symmetric links, [18] proposed two slave clocks, triggering the clock counters in the slaves. Figure 2 (recalled from [18] and slightly modified) describes the PTP messaging required by the dual slave clocks method, where is the time when the th Sync message is sent. and are the times that slave clocks 1 and 2 measure at the receive of the th Sync message, respectively. Notably, the dual slave clocks method requires only one-way time transfer. Clock 2 runs at a normalized clock period of which is normalized to that of the master clock, where is the normalized clock frequency difference between slave clock 2 and the master clock (i.e., normalized clock skew). A D-type flip-flop is employed to divide the frequency by two. Therefore, slave clock 2 runs at a normalized clock period of , which is half of the slave clock 1 running at a normalized clock period of . Both slave clocks have the same offset, which is assumed to be . Thus, according to [18] we may express the times , as

where

(i) and are random jitters of slave clock 1 and 2 within the th period, respectively.(ii)The random jitter is assumed to be an independent Gaussian random variable (RV) with zero mean and variance .(iii)Since slave clock 2 is derived from slave clock 1 and the period of slave clock 2 is twice that of slave clock 1, the variance of is due to accumulation of random jitters at two positive edges of slave clock 1.

Here ends the brief presentation of the dual slave clock method in a slave as presented in [18].

Figure 2: Synchronization of dual slave clocks.

The uplink delay can be represented as a constant value with an additive zero mean random noise called PDV. The PDV as already was mentioned results from the behaviour of switches or routers in the packet network [24]. In the following, the constant value is the we saw in [18], and will be denoted by the PDV of the th cycle. Therefore, (8) must be written in the following form:

where is assumed to be fractional Gaussian noise (fGn) with zero mean with Hurst exponent and variance .

3. ML Clock Skew Estimation for Fractional Gaussian Noise PDV

ML estimators have a number of attractive characteristics. For example, they are asymptotically unbiased and efficient.

In this section we derive the ML estimation and the CRLB for the clock skew based on the dual slaves clock model presented in [18] where the PDV is modeled with fGn unlike how it was done in [18].

In Section 3.1 we derive the ML estimation for the clock skew under the assumption of fGn PDV with Hurst exponent in the range of .

In Section 3.2 we derive an expression for the Cramer-Rao (CR) lower bound on the variance of the clock skew.

3.1. ML Estimation of Clock Skew

Theorem 1. For the following assumptions

(i)the Hurst exponent is in the range of ,(ii),

is assumed to be fGn with zero mean with Hurst exponent and variance .

The ML estimation for the clock skew may be defined as

where is defined by

denotes Gamma function and is the number of observation.

Comments

(i) denotes the difference between two consecutive samples.(ii)The system equations are (9) and (10) instead of (8) in order to describe the system more closely to a real system with PDV.

Proof. Subtracting (10) from (9) yields

Subtracting (13) of the th Sync message from that of the th Sync message indicates that

which can be also written as

Thus, may be expressed as

where according to [42]

According to [49, 50], the probability density function (PDF) of vector (of length ) may be given by

where denotes the transpose of and denotes the determinant of the matrix .
is the covariance matrix of the vector :

The elements of are given by [50, 51]

where, according to [51, 52],

Note that so .
can be decomposed as was done in [49, 50]

where

According to [50] we can write

where

Next, we recall from [50] the expression for the log-likelihood function:

The derivative of the log-likelihood function (26) with respect to yields

Note that [50] has taken the derivative of (26) with respect to . Thus, our final expression is very different.
From (24) and (25) we can write

By substituting (28) into (27) we obtain

According to [50], the expression can be written as

where

Then, we substitute (30) into (29) and obtain

Let us define for simplicity the next two variables:

Inserting (34) and (35) into (33) gives

Next, the derivative of the log-likelihood function (36) with respect to is carried out and set to zero

Appendix A shows in detail the derivation of the log-likelihood function (36) and the calculation of (37). In Appendix C we give the sum limitation of (36).
According to Appendices A and C we obtain the ML estimation of the clock skew as written in (11), where is obtained in Appendix B as given in (12).
This completes our proof.

3.2. CRLB for the ML Estimation of Clock Skew

Theorem 2. For the following assumption
(i)   , where is the fixed interval between consecutive PTP Sync messages.
The CRLB of the clock skew ML estimation may be defined as

where denotes the variance of and is defined by

and is defined by

Proof. The inequality of the Cramer-Rao bound is given by

First, we will calculate the second derivative of the log-likelihood function.
From (36) and Appendix C we can write the derivative of the log-likelihood function with respect to as

Inserting (16) into (42) gives

Therefore, according to Appendix D, the second derivative of the log-likelihood function with respect to may be written as

where the full proof is provided in Appendix D.
Next, we will derive the mean of (44) and according to Appendix D we may write

Finally, by substituting (45) into (41) and using some properties obtained in Appendix D, we come back to the expression for the CRLB as given in (38).
This completes our proof.

4. Simulation

In this section we show the MSE performance comparison between our proposed expression for the clock skew (11) compared with the conventional dual slaves clock method proposed in [18]. We also test the performance of our new proposed expression for the clock skew (11) for two different values of and . In our simulation the Sync interval and the clock skew have the same values as in [18]: the Sync interval  sec and the clock skew  ppm.

Figures 3 and 4 show the MSE performance comparison between the proposed approach (11) and the conventional dual slaves clock method given in equation in [18] for and , respectively. Namely, the MSE performance comparisons (Figures 3 and 4) are carried out for the case where the PDV is modeled with fGn. The standard deviation is the same value as was taken in [18] for the random jitter ( ps). Figures 3 and 4 show that our proposed expression (11) for the clock skew outperforms (in the MSE point of view) the dual slaves clock method given in equation in [18] for both cases and . Thus we may conclude that the proposed algorithm for the clock skew [18] may not fit for a realistic Ethernet network where PDV exists. But it should be pointed out that the proposed algorithm for the clock skew [18] is very simple compared to our proposed approach (11). As a matter of fact the computational complexity of (11) is while the computational complexity of [18] is only (where is defined as and “” is a constant).

Figure 3: MSE comparison between the proposed ML estimator for the clock skew (11) to the conventional dual slave method ( in [18]) as a function of the observation number . The standard deviation is  ps.

Figure 4: MSE comparison between the proposed ML estimator for the clock skew (11) to the conventional dual slave method ( in [18]) as a function of the observation number . The standard deviation is  ps.

Figures 5–8 show the MSE performance obtained from our proposed clock skew estimator (11) as a function of for two different values of and .

Figure 5: MSE and CRLB of the clock skew as a function of the observation number . The standard deviation of the PDV is  ps.

Figure 6: MSE and CRLB of the clock skew as a function of the observation number . The standard deviation of the PDV is  ps.

Figure 7: MSE and CRLB of the clock skew as a function of the observation number . The standard deviation of the PDV is  ps.

Figure 8: MSE and CRLB of the clock skew as a function of the observation number . The standard deviation of the PDV is  ps.

Since the simulation tool (MATLAB) can not calculate the Gamma function for an argument higher than 175, we can not run the simulation for . Therefore, we have to reduce the value for in order to be able to show results also for higher values of the Hurst exponent for which the MSE approaches the CRLB. Nevertheless, in order to show it, Figures 5 and 6 were produced with a lower value for (than chosen in Figures 7-8) in order to show that the MSE approaches the CRLB also for higher values of the Hurst exponent as expected. According to Figures 5-6 we can see that the MSE declines and approaches the CRLB (38) as increases. However, according to Figures 7 and 8 we can see that even for the higher values for the MSE declines and goes in the right direction towards the CRLB, which means that for greater values for we will see a complete match between the MSE to the CRLB.

5. Conclusion

In this paper we developed a new clock skew method for the case where the PDV is modeled as fGn, as proposed in the literature for a realistic Ethernet network. Simulation results indicate that our new proposed method for the clock skew leads the system with a lower MSE compared with the literature proposed dual slave clocks in a slave. In this paper we derived also the CRLB ( dependent) for our proposed clock skew method and have shown that the simulated MSE obtained from our proposed method achieves asymptotically our new derived CRLB.

Appendices

A. Log-Likelihood Function Derivation

In this appendix we derive in detail the derivative of the log-likelihood function with respect to .

We also show the calculation of the log-likelihood function maximum value.

Equation (36) determines

Since is not dependent on , we can write

Then, we substitute (16) into (A.2) and obtain

Now, we will derive the internal expression with respect to and get

Next, we set (A.4) to zero in order to find the log-likelihood function maximum value:

After some basic algebraic operations we get

Hence, by substituting (17) into (A.6) we obtain

Since we can say that

Now, we will multiply the numerator and denominator in the first term of (A.7) by :

Inserting (A.8) into (A.9) gives

Since the random jitters and are assumed to be independent we will write

Inserting (A.11) into (A.10) gives

Since the random jitter is assumed to be an independent Gaussian RV with zero mean and variance of we can write

Therefore

where only when . because occurs times. because or occurs times.

We will use the same technique regarding :

Therefore

where only when . because occurs times. because or occurs times.

Inserting (A.14) and (A.16) into (A.12) gives

Finally, according to Appendix C we may write (A.17) with the proper sum limitation:

B. Derivation

In this appendix we derive in detail the expression for .

Earlier we defined in (34):

By setting (32) into (B.1) we obtain

Now, by substituting (21) into (B.2) and performing some algebraic operation we obtain

We will write with the Gamma function instead of the factorial function because .

The main character of the Gamma function is

Therefore

Finally, we will arrange the expression for more elegantly

This is the final form of mentioned in (12).

C. Log-Likelihood Function First Derivative Sum Limitation

In this Appendix we derive in detail the sum limitations: .

We can get the above limitation by two different ways:

(i)examining the limitations on the Gamma function arguments that exist in ,(ii)examining the limitation on the expression mentioned in [50] that built up the expression.

The first way:

(i) The Gamma function has simple poles at .

Now, under the above knowledge, we examine all the Gamma function arguments existing in (B.6):

(i.A) always exists;(i.B)   exists under (1.J);(i.C) always exists;(i.D) exists under (1.K);(i.E) always exists;(i.F) always exists;(i.G) always exists;(i.H) always exists;(i.I) always exists;(i.J);(i.K);(i.L) exists under (1.J);(i.M) exists under (1.K);(i.N) always exists;(i.O) always exists.

From (i.J) and (i.K) we get the limitation: .

Note that (i.A)–(i.I) and (i.L)–(i.O) are always true for the abovementioned limitations.

Therefore, the conclusion from the first way is

The second way:

(ii) From the limitation mentioned in [50, page 2981], we can learn that is only defined for and .

Therefore

(ii.A)   and ,(ii.B) and ,(ii.C) and ,(ii.D) and .

From (ii.A) and (ii.B) and from the fact that we get that .

From (ii.C) and (ii.D) we get the limitation: .

In conclusion we get the same limitations as in (C.1).

D. Log-Likelihood Function-Second Derivative

In this Appendix we derive in detail the Cramer-Rao lower bound (CRLB).

The inequality of the Cramer-Rao bound as given in (41) is

In order to calculate the second derivative of the log-likelihood function we will write the first derivative of the log-likelihood function with respect to .

We can write it from (36)

Then, we substitute (16) for the th sample and (16) for the th sample into (D.2) and obtain

Next, we will take the derivative according to the product role and get

Now, let us calculate the second derivative of the log-likelihood function by deriving (D.4) with respect to

Inserting (16) into (D.5) we obtain

Therefore

Now, we will calculate the mean of the second derivative of the log-likelihood function:

Because and are deterministic values we can say

Since and are assumed to be with zero mean, we can write

Now, by substituting (D.10) into (D.9) we obtain

According to Appendices E and F and by substituting (20) into (D.11) we obtain

where

Then, by substituting (35) and (D.13) into (D.1) we obtain

Let us divide the numerator and denominator by

In order to find we will set into (23) and obtain

Similarly, in order to find we will set into (20) and obtain

By substituting (D.19) into (D.18) we obtain

Then, we substitute (D.20) into (D.17) and obtain

Let us divide the numerator and denominator by

Finally, by substituting (D.14) into (D.22) we obtain

According to Appendices E and G and the above equation (D.23) the CRLB is defined by

where

and is obtained by

E. Calculation

In this Appendix we show in detail how to calculate .

According to equation given in [50], the covariance function of a fractionally differenced Gaussian noise process (fdGn) is given by

Therefore, we can write

On the other hand, we defined in this paper (assumption in subsection “ML Estimation of Clock Skew”) , as

By substituting (E.3) into we can say that

Next, by substituting (E.1) into (E.4) we obtain

By comparing (E.2) with (E.5) we get

Finally, by substituting (D.14) into (E.6), we get

Therefore, the definition of is given by

where

F. Calculation

In this Appendix we show in detail how to calculate .

We will calculate according to the Rescaled-Range () analysis.

According to [51, 52] in order to estimate the Hurst exponent of a time series ( denotes vector of data set), we should execute the following process.

First, the time series of full length is divided into a number of shorter time series of length .

The average Rescaled-Range is then calculated for each value of .

The Rescaled-Range is calculated as follows.

(i) Mean calculation:

The mean of the time series is calculated as follows:

(ii) Adjusted mean series calculation:

(iii) Cumulative deviate series calculation:

(iv) Range series calculation:

(v) Standard deviation series calculation:

(vi) Average over all the partial time series of length of the Rescaled-Range :

And finally are (vii) plotting the logarithm as a function of and fitting a straight line according to the linear regression method. The slope of the line gives .

Figure 9 shows the plot of the Rescaled-Range analysis as a function of (where is defined as the log function with base of 2) for the case of .

Figure 9: The analysis in case of Hurst exponent . The obtained zeta Hurst exponent is .

Operating the above analysis on several data sets (corresponding to several values of ) and by taking the mean value gives the following results:

The values for are denoted by and the values for are denoted by .

By a second-order polynominal fit preformed in Excel, we find the following strong () connection between and :

where is the Pearson product-moment correlation coefficient and the value can be interpreted as the proportion of the variance in attributable to the variance in .

is a measure of the dependence between two variables, giving a value between and −1 inclusive, where 1 is total positive correlation. It is widely used in the sciences as a measure of the degree of dependence between two variables.

can be calculated by the following equation:

where

(i) and are the mean of the and data sets, respectively.(ii) and are the standard deviation of the and data sets, respectively.

G. Derivation

In this Appendix we derive in detail the definition of and its sum limitations: .

Inserting (E.8) into (D.23) gives

Therefore, is defined by

where is defined as

Next we will explain the sum limitations of .

First, the sum limitations contains at least the same limitations as in (shown in Appendix C) because is a product of with additional term. Therefore, the first limitation will be

Second, we will examine the limitations on the Gamma function arguments existing in the second term of :

As already mentioned in this paper, the Gamma function has simple poles at ; therefore

(A) always exists;(B) always exists under ;(C) always exists;(D) always exists under ;(E) always exists;(F) always exists;(G) always exists;(H) always exists under ;(I) always exists under ;(J) always exists under ;(K) always exists under ;(L) always exists under ;(M) always exists under ;

Finally, we will examine the terms in the (G.5) denominator:

(i) always exists; .(ii) always exists; .

Therefore, the conclusion is that has the same limitations like .

Conflict of Interests

The authors declare that there is no conflict of interests regarding the publication of this paper.

Acknowledgment

The authors would like to thank the anonymous reviewers for their helpful comments.

References

O. Pallars, J. Cadena-Muoz, N. Carreras, D. M. Toma, and J. del-Ro-Fernndez, “Contribution to COBS synchronization with PTP IEEE std. 1588,” in Proceedings of the MTS/IEEE OCEANS—Bergen, pp. 1–4, Bergen, Norway, June 2013.

B. Ferencz and B. Kovcshzy, “Hardware assisted COTS IEEE 1588 solution for x86 linux and its performance evaluation,” in Proceedings of the International IEEE Symposium on Precision Clock Synchronization for Measurement Control and Communication (ISPCS '13), pp. 47–52, Lemgo, Germany, September 2013.

P. Ferrari, A. Flammini, S. Rinaldi, and G. Prytz, “High availability IEEE 1588 nodes over IEEE 802.1aq shortest path bridging networks,” in Proceedings of the International IEEE Symposium on Precision Clock Synchronization for Measurement Control and Communication (ISPCS '13), pp. 35–40, 2013.

N. Moreira, A. Astarloa, J. Lazaro, A. Garcia, and E. Ormaetxea, “IEEE 1588 Transparent Clock architecture for FPGA-based network devices,” in Proceedings of the IEEE 22nd International Symposium on Industrial Electronics (ISIE '13), pp. 1–6, IEEE, Taipei, Taiwan, May 2013.

J. Jeon, D. Kim, and D. Lee, “IEEE1588-based clock synchronization for embedded networked system with sRIO,” in Proceedings of the International Conference on Information and Communication Technology Convergence (ICTC '13), pp. 843–845, IEEE, Jeju, Republic of Korea, October 2013.

Q. Yang, D. An, and W. Yu, “On time desynchronization attack against IEEE 1588 protocol in power grid systems,” in Proceedings of the IEEE EnergyTech 2013 Conference, pp. 1–5, May 2013.

H. Flatt, S. Schriegel, and J. Jasperneite, “Reliable synchronization accuracy in ieee 1588 networks using device qualification with standard test patterns,” in Proceedings of the 7th International IEEE Symposium on Precision Clock Synchronization for Measurement, Control and Communication (ISPCS '13), pp. 107–112, IEEE, Lemgo, Germany, September 2013.

A. Mahmood and R. Exel, “Servo design for improved performance in software timestamping-assisted WLAN synchronization using IEEE 1588,” in Proceedings of the IEEE 18th International Conference on Emerging Technologies and Factory Automation (ETFA '13), pp. 1–8, September 2013.

N. Moreira, A. Astarloa, and U. Kretzschmar, “SHA-3 based Message Authentication Codes to secure IEEE 1588 synchronization systems,” in Proceedings of the 39th Annual Conference of the IEEE Industrial Electronics Society (IECON '13), pp. 2323–2328, November 2013.

J. del Río, D. M. Toma, S. Shariat-Panahi, Y. Auffret, and A. Mánuel-Lázaro, “Smart IEEE-1588 GPS clock emulator for cabled ocean sensors,” IEEE Journal of Oceanic Engineering, vol. 39, no. 2, pp. 269–275, 2014.

J. Peng, L. Zhang, and D. McLerno, “On the clock offset estimation in an improved IEEE 1588 synchronization scheme,” in Proceedings of the 19th European Wireless Conference (EW '13), pp. 1–6, Guildford, UK, April 2013.

ITU-T Recommendation, “G.8262 Timing Characteristics of Synchronous Ethernet Equipment Slave Clock (EEC),” 2007

D. L. Mills, “Precision synchronization of computer network clocks,” ACM SIGCOMM Computer Communication Review, vol. 24, no. 2, pp. 28–43, 1994.

IEEE Instrumentation and Measurement Society, “IEEE 1588 standard for a precision clock synchronization protocol for networked measurement and control systems IEEE Std,” Tech. Rep. 1588-2008 version 2, 2008.

P. Yang, H. Tang, X. Chen, and M. Tian, “New algorithm for IEEE 1588 time synchronization under the presence of significant delay variation,” in Proceedings of the 3rd International Conference on Consumer Electronics, Communications and Networks (CECNet '13), pp. 419–422, November 2013.

Request for Comment 5905 (RFC), Network Time Protocol Version 4: Protocol and Algorithms Specification, Internet Engineering Task Force (IETF), 2010.

J. Ridoux and D. Veitch, “Ten microseconds over LAN, for free (extended),” IEEE Transactions on Instrumentation and Measurement, vol. 58, no. 6, pp. 1841–1848, 2009.

W.-L. Chin and S.-G. Chen, “IEEE 1588 clock synchronization using dual slave clocks in a slave,” IEEE Communications Letters, vol. 13, no. 6, pp. 456–458, 2009.

P. Ferrari, A. Flammini, D. Marioli, and A. Taroni, “IEEE 1588-based synchronization system for a displacement sensor network,” IEEE Transactions on Instrumentation and Measurement, vol. 57, no. 2, pp. 254–260, 2008.

T. Cooklev, J. C. Eidson, and A. Pakdaman, “An implementation of IEEE 1588 over IEEE 802.11b for synchronization of wireless local area network nodes,” IEEE Transactions on Instrumentation and Measurement, vol. 56, no. 5, pp. 1632–1639, 2007.

A. Vallat and D. Schneuwly, “Clock synchronization in telecommunications via PTP (IEEE 1588),” in Proceedings of the IEEE International Frequency Control Symposium Joint with the 21st European Frequency and Time Forum (FCS '07), pp. 334–341, June 2007.

J. Eidson, IEEE-1588 Standard for a Precision Clock Synchronization Protocol for Networked Measurement and Control Systems—A Tutorial, Agilent Technologies, Santa Clara, Calif, USA, 2005.

P. Ferrari, A. Flammini, D. Marioli, and A. Taroni, “IEEE 1588-based synchronization system for a displacement sensor network,” IEEE Transactions on Instrumentation and Measurement, vol. 57, no. 2, pp. 254–260, 2008.

ITU-T G.8261/Y.1361, Timing and synchronization aspects in packet networks, 2006

M. Li, W. Jia, and W. Zhao, Simulation of Long-Range Dependent Traffic and a Simulator of TCP Arrival Traffic, World Scientific Publishing, River Edge, NJ, USA, 2001.

M. Li, C.-H. Chi, and D. Long, “Fractional Gaussian noise: a tool of characterizing traffic for detection purpose,” in Content Computing, vol. 3309 of Lecture Notes in Computer Science, pp. 94–103, Springer, Berlin, Germany, 2004.

O. Cappe, E. Moulines, J.-C. Pesquet, A. P. Petropulu, and X. Yang, “Long-range dependence and heavy-tail modeling for teletraffic data,” IEEE Signal Processing Magazine, vol. 19, no. 3, pp. 14–27, 2002.

M. T. Lucas, D. E. Wrege, B. J. Dempsey, and A. C. Weaver, “Statistical characterization of wide-area IP traffic,” in Proceedings of the 6th International Conference on Computer Communications and Networks, vol. 1095, pp. 442–447, IEEE, Las Vegas, Nev, USA, September 1997.

J. Jusak and R. J. Harris, “Study of UDP-based Internet traffic: long-range dependence characteristics,” in Proceedings of the Australasian Telecommunication Networks And Applications Conference (ATNAC '11), pp. 1–7, November 2011.

M. Gospodinov and E. Gospodinova, “The graphical methods for estimating Hurst parameter of self-similar network traffic,” in Proceedings of the International Conference on Computer Systems and Technologies (CompSysTech '05), Varna, Bulgaria, June 2005.

M. Li, “Fractional Gaussian noise and network traffic modeling,” in Proceedings of the 8th WSEAS International Conference on Applied Computer and Applied Computational Science, 2009, http://www.wseas.us/e-library/conferences/2009/hangzhou/ACACOS/ACACOS03.pdf.

M. Li, S. C. Lim, and H. Feng, “Generating traffic time series based on generalized cauchy process,” in Computational Science—ICCS 2007: 7th International Conference, Beijing, China, May 27–30, 2007, Proceedings, Part I, vol. 4487 of Lecture Notes in Computer Science, pp. 374–381, Springer, Berlin, Germany, 2007.

B. Tsybakov and N. D. Georganas, “Self-similar processes in communications networks,” IEEE Transactions on Information Theory, vol. 44, no. 5, pp. 1713–1725, 1998.

M. Li, W. Zhao, W. Jia, D. Long, and C.-H. Chi, “Modeling autocorrelation functions of self-similar teletraffic in communication networks based on optimal approximation in Hilbert space,” Applied Mathematical Modelling, vol. 27, no. 3, pp. 155–168, 2003.

V. Paxson and S. Floyd, “Wide area traffic: the failure of poisson modeling,” IEEE/ACM Transactions on Networking, vol. 3, no. 3, pp. 226–244, 1995.

V. Paxon, “Approximate synthesis of fractional Gaussian noise for generating self-similar network traffic,” Computer Communications Review, vol. 27, no. 5, pp. 5–18, 1997.

H.-D. J. Jeong, J.-S. R. Lee, D. McNickle, and K. Pawlikowski, “Distributed steady-state simulation of telecommunication networks with self-similar teletraffic,” Simulation Modelling Practice and Theory, vol. 13, no. 3, pp. 233–256, 2005.

S. Ledesma and D. Liu, “Synthesis of fractional gaussian noise using linear approximation for generating self-similar network traffic,” Computer Communication Review, vol. 30, pp. 4–17, 2000.

M. W. Garrett and W. Willinger, “Analysis, modeling and generation of self-similar VBR video traffic,” in Proceedings of the Conference on Communications Architectures, Protocols and Applications (SIGCOMM '94), pp. 269–280, London, UK, September 1994.

M. Li and C.-H. Chi, “A correlation-based computational model for synthesizing long-range dependent data,” Journal of the Franklin Institute, vol. 340, no. 6-7, pp. 503–514, 2003.

M. Li, “An approach to reliably identifying signs of DDOS flood attacks based on LRD traffic pattern recognition,” Computers and Security, vol. 23, no. 7, pp. 549–558, 2004.

K.-L. Noh, Q. M. Chaudhari, E. Serpedin, and B. W. Suter, “Novel clock phase offset and skew estimation using two-way timing message exchanges for wireless sensor networks,” IEEE Transactions on Communications, vol. 55, no. 4, pp. 766–777, 2007.

T. Trump, “Estimation of clock skew in telephony over packet switched networks,” in Proceedings of the IEEE Interntional Conference on Acoustics, Speech, and Signal Processing, pp. 2605–2608, June 2000.

X. Cao, F. Yang, X. Gan et al., “Joint estimation of clock skew and offset in pairwise broadcast synchronization mechanism,” IEEE Transactions on Communications, vol. 61, no. 6, pp. 2508–2521, 2013.

M. Leng and Y.-C. Wu, “Low-complexity maximum-likelihood estimator for clock synchronization of wireless sensor nodes under exponential delays,” IEEE Transactions on Signal Processing, vol. 59, no. 10, pp. 4860–4870, 2011.

M. Leng and Y.-C. Wu, “On clock synchronization algorithms for wireless sensor networks under unknown delay,” IEEE Transactions on Vehicular Technology, vol. 59, no. 1, pp. 182–190, 2010.

M. Leng and Y.-C. Wu, “On joint synchronization of clock offset and skew for Wireless Sensor Networks under exponential delay,” in Proceedings of IEEE International Symposium on Circuits and Systems (ISCAS '10), pp. 461–464, Paris, France, May-June 2010.

Q. M. Chaudhari, E. Serpedin, and K. Qaraqe, “On maximum likelihood estimation of clock offset and skew in networks with exponential delays,” IEEE Transactions on Signal Processing, vol. 56, no. 4, pp. 1685–1697, 2008.

T. Lundahl, W. J. Ohley, S. M. Kay, and R. Siffert, “Fractional Brownian motion: a maximum likelihood estimator and its application to image texture,” IEEE Transactions on Medical Imaging, vol. 5, no. 3, pp. 152–161, 1986.

M. Deriche and A. H. Tewfik, “Maximum likelihood estimation of the parameters of discrete fractionally differenced Gaussian noise process,” IEEE Transactions on Signal Processing, vol. 41, no. 10, pp. 2977–2989, 1993.

J. R. M. Hosking, “Fractional differencing,” Biometrika, vol. 68, no. 1, pp. 165–176, 1981.

T. Dieker, Simulation of Fractional Brownian Motion, 2015.