Maximum-Likelihood Estimator of Clock Offset between Nanomachines in Bionanosensor Networks

Recent advances in nanotechnology, electronic technology and biology have enabled the development of bio-inspired nanoscale sensors. The cooperation among the bionanosensors in a network is envisioned to perform complex tasks. Clock synchronization is essential to establish diffusion-based distributed cooperation in the bionanosensor networks. This paper proposes a maximum-likelihood estimator of the clock offset for the clock synchronization among molecular bionanosensors. The unique properties of diffusion-based molecular communication are described. Based on the inverse Gaussian distribution of the molecular propagation delay, a two-way message exchange mechanism for clock synchronization is proposed. The maximum-likelihood estimator of the clock offset is derived. The convergence and the bias of the estimator are analyzed. The simulation results show that the proposed estimator is effective for the offset compensation required for clock synchronization. This work paves the way for the cooperation of nanomachines in diffusion-based bionanosensor networks.


Introduction
Recent advances in the fields of nanotechnology, electronic technology and biology have enabled the development of bionano-devices manufactured in a scale ranging from one to a hundred nanometers [1]. At this scale, a nano-device, also called a bionanosensor or nanomachine, is considered to be the most basic functional device, able to perform tasks at the nano-level, such as computing, data storage, sensing or actuation. Like ordinary sensors and wireless sensor networks, nanomachines can also be interconnected, forming a bionanosensor network [2], to execute more complex tasks and to expand their range of operations [3]. Bionanosensor networks have a broad range of potential applications. For example, in the area of biomedical engineering, bionanosensor networks can be applied in intra-body health monitoring, drug delivery systems, etc. In industrial and consumer goods applications, bionano-sensor networks can be used for the development of intelligent functionalized materials, new manufacturing processes, distributed quality control procedures, etc. [4].
Molecular communication, a kind of bio-inspired communication, is a very popular and promising information exchange technique in the recent nano-communication and bionanosensor network studies [5][6][7][8]. A typical molecular communication system is shown in Figure 1. A transmitter nanomachine sends messages by releasing the information molecules. The information molecules propagate directionally or diffuse randomly in the environment, and reach a receiver nanomachine. The receiver nanomachine then decodes the information molecules into the original message. Since the molecular communication paradigm is different from the conventional radio communication at the physical layer such as channel and signal propagation, the corresponding communication techniques and networking protocols need to be reconsidered and redesigned.
information exchange technique in the recent nano-communication and bionanosensor network studies [5][6][7][8]. A typical molecular communication system is shown in Figure 1. A transmitter nanomachine sends messages by releasing the information molecules. The information molecules propagate directionally or diffuse randomly in the environment, and reach a receiver nanomachine. The receiver nanomachine then decodes the information molecules into the original message. Since the molecular communication paradigm is different from the conventional radio communication at the physical layer such as channel and signal propagation, the corresponding communication techniques and networking protocols need to be reconsidered and redesigned.

Figure 1. A typical molecular communication system in a bionnanosensor network.
Clock synchronization is a very important issue for all kinds of communication networks, including nano-communication and nanosensor networks. It is the foundation of the distributed collaboration among a set of cooperative nanomachines. For example, in an immune defense bionnanosensor network, the nanomachines are envisioned to release antibody molecules at the same time to effectively kill the pathogen, or in a nanosensor network-based monitoring scenario, clock synchronization is the basis for the interpretation of the sensed data from different nanomachines. Clock synchronization is a very important issue for all kinds of communication networks, including nano-communication and nanosensor networks. It is the foundation of the distributed collaboration among a set of cooperative nanomachines. For example, in an immune defense bionnanosensor network, the nanomachines are envisioned to release antibody molecules at the same time to effectively kill the pathogen, or in a nanosensor network-based monitoring scenario, clock synchronization is the basis for the interpretation of the sensed data from different nanomachines. However, the slowness and randomness of the propagating molecules represent a big challenge to clock synchronization in bionanosensor networks. The relevant studies of the clock synchronization issues in the wireless sensor networks field can be the good references for the study of clock synchronization in bionanosensor networks [9,10], while those solutions cannot be directly applied to solving the clock synchronization problems in the bionnanosensor networks. The solutions in [11,12] proposed an averaging algorithm for the clock offset compensation for the synchronization in the wireless sensor networks. However, the design ignores the propagation delay of the information exchange between the sensor nodes. Since in bionnanosensor networks, the information molecules propagate much more slowly than electromagnetic waves, the propagation delay is much larger, which cannot be ignored in any algorithm design for clock synchronization. The timing-sync protocol for sensor networks (TPSN) [13] and the recursive time synchronization protocol (RTSP) [14] for wireless sensor networks considered the propagation delay in the design, however, they assume that the propagation delays are bi-directionally symmetrical, while in molecular communications, the molecules move with strong randomness, so the assumption of symmetrical propagation delays is not appropriate for molecular communications. References [15,16] modeled the random propagation delay following a Gaussian distribution. References [16][17][18][19][20] modeled the random delay as an exponential delay and a bivariate exponential delay, but in molecular communications, the propagation delay is normally assumed to follow an inverse Gaussian distribution [21]. In [22], an approach for clock offset estimation has been proposed which is robust to the distribution of the network delays, however, it has high computational complexity and no closed-form expression of the estimator can be obtained.
The oscillation and synchronization of bionanosensor networks have been investigated in recent years. In [23,24], the authors used a bacterial quorum sensing mechanism to synchronize nano-machines. Molecules, called inducers, can be released by one nanomachine, and trigger another nanomachine to release the same self-induced molecules. When the concentration of the inducer particles in the environment reaches a certain threshold at a moment in time, the entire nanonetwork can achieve clock synchronization. In [25,26], clock synchronization was realized by inhibitory molecules. These molecules are released by one transmitter, and inhibit the release of the similar molecules from another transceiver in the nanonetwork. When the concentration of molecules falls below a certain threshold, those molecules can be released again. The pattern of the release pulses of the inhibitory molecules provides the synchronization. However, the abovementioned works try to synchronize the oscillation period rather than align the time. It is crucial and indispensable for the nanomachines to synchronize the time in various nanonetwork applications. In [27], the authors proposed a blind algorithm for such synchronization using non-decision directed maximum likelihood. The clock sequence is calculated by the receiver based on the analysis of the molecular channel delay. However, the authors designed receiver sampling sequences without any further discussion on how to synchronize the time between the transmitter and the receiver.
In this paper, based on the statistical delay model of the molecular diffusion, a two-way message exchange model for clock synchronization is proposed. A maximum likelihood estimator of the clock offset between two nanomachines is designed. We assume that the nanomachines have their own time and the aim of this work is to estimate the time difference between nanomachines. This is the major difference of this paper from the other literature [6][7][8][9][10]. To the best of our knowledge, this is the first work focusing on estimating the clock offsets for diffusion-based nano-communication. The simulation results validate the effectiveness of the proposed algorithm. This work is an important step and the foundation for the distributed cooperation of nanomachines.
The rest of the paper is organized as follows: Section 2 introduces the system model. The maximum likelihood estimator and the analytical analysis are given in Section 3. The Cramer Rao Lower Bound is derived in Section 4. Section 5 presents the simulation results and discussions. Finally, Section 6 concludes the paper.

System Model
In bionanosensor networks, it is assumed that every individual nanomachine has its own clock. The clock in the nanomachine is analogous to the clock in electronic devices [28]. A clock in an electronic device is considered as a measurement device consisting of an oscillator and memory. Accordingly, the nanomachine components are designed similarly. The memory and processing functionalities in nanomachines are mentioned and designed in [4,29,30]. For example, the nanomachine can use the increase or decrease Ca 2+ levels, high or low levels of the pH value, or some three dimensional structure to implement the memory function. There are several ways to implement the oscillator. A circadian oscillator is proposed in [31]. In references [32,33] translation and transcription in gene expression is used to generate the oscillation. A molecular phase locked loop (PLL) is proposed in [28] which can be used for the oscillator in a nanomachine. The clock reading of node i at time t, C i ptq, is defined as in Equation (1), where β i is the clock offset. In the ideal situation, β i is equal to 0, and C i ptq is equal to the real time t: In order to achieve a distributed cooperation among nanomachines, clock synchronization is required. The objective of clock synchronization is to estimate and compensate the clock offset. For this purpose, one nanomachine can send its clock readings to another nanomachine. Since the message exchange in the nanonetwork is through molecular communication, the information molecules need a quite long period of time to propagate from the transmitter to the receiver. Therefore the propagation delay estimation is necessary.
In a molecular communication scenario, we consider a clock synchronization process between two nanomachines, node A and node B. Assume that node A and node B have independent times.
The clock offset between them is supposed to be β 0 , where β 0 " β B´β A . Node A is synchronized with node B by a two-way timing message exchange mechanism shown in Figure 2. Take the ith round of message exchange as an example. Node A sends a clock synchronization request message. At the same time it records the timing on its own clock as T 1,i . Node B records its time T 2,i at the reception of the message, and replies to node A at T 3,i . The reply message contains the time stamps T 2,i and T 3,i . Then, node A records the reception time of node B's reply as T 4,i . T 1,i and T 4,i are the time stamps recorded by the clock of node A, while T 2,i and T 3,i are those recorded by node B. After n rounds of message exchange, node A obtains a set of time stamps Node A then tries to use these time stamps to estimate the propagation delay of the information molecules. The above procedure can be modeled as: where β 0 denotes the clock offset of node B with respect to node A. X i and Y i are the random variables of the propagation delays from node A to node B and from node B to node A, respectively.
Sensors 2015, 15 5 where 0 denotes the clock offset of node B with respect to node A. and are the random variables of the propagation delays from node A to node B and from node B to node A, respectively. Molecular diffusion is characterized by Brownian motion. It is assumed in this paper that the information is modulated based on the type of the molecules. Each molecule can carry multiple information bits [34] so that the time stamp message can be encoded into an individual molecule. Normally, a Gaussian distribution is used to estimate the distance of the Brownian motion in a fixed period. The inverse Gaussian distribution is used to describe the distribution of the duration a Brownian motion takes to reach a fixed positive distance [35]. The probability density function of an inverse Gaussian distribution random variable can be expressed as: where t is the random propagation delay, known as the first arrival time [21]. is the mean of t. is the shape parameter. According to [36], = and =  Molecular diffusion is characterized by Brownian motion. It is assumed in this paper that the information is modulated based on the type of the molecules. Each molecule can carry multiple information bits [34] so that the time stamp message can be encoded into an individual molecule. Normally, a Gaussian distribution is used to estimate the distance of the Brownian motion in a fixed period. The inverse Gaussian distribution is used to describe the distribution of the duration a Brownian motion takes to reach a fixed positive distance [35]. The probability density function of an inverse Gaussian distribution random variable can be expressed as: where t is the random propagation delay, known as the first arrival time [21]. µ is the mean of t. λ is the shape parameter. According to [36], µ " d v and λ " d 2 2D , where d represents the distance between the nanomachines, v is the medium velocity, and D represents the diffusion coefficient. X i and Y i are modeled as independent and identical distributed (i.i.d.) inverse Gaussian random variables with a common shape parameter λ. The goal is to estimate the clock offset β 0 based on the observation of a set of time stamps.

Maximum-Likelihood Estimator
In this section, we propose an approach using maximum-likelihood estimation (MLE) to estimate the clock offset β 0 . To derive the MLE, we rewrite Equations (2) and (3) as Equations (5) and (6): Since tX i , Y i u N i"1 are i.i.d. inverse Gaussian random variables, denoting their probability density function as tp pX i q , p pY i qu N i"1 , the likelihood function of T 1, i , T 2, i , T 3, i , T 4, i ( N i"1 is given by ś N i"1 p px i q p py i q. x i and y i are the outcomes of X i and Y i . Equations (5) and (6) into the likelihood function, we have the expression: In the system model, µ and λ are unknown by the nanomachines. They are considered as the nuisance parameters. To estimate β 0 , µ and λ should be replaced by the expression of T 1,i , T 2,i , T 3,i , T 4,i ( N i"1 and β 0 . A new random variable, X i`Yi , is introduced. By this way β 0 can be eliminated. Since X i and Y i are i.i.d. inverse Gaussian random variables, X i`Yi " T 4,i´T3,iT 2,i´T1,i also follows the inverse Gaussian distribution, with a mean of 2µ and a shaping parameter of 4λ [35]. The likelihood function can be expressed as: The MLE of µ and λ can be obtained by differentiating the logarithm of Equation (8) with respect to µ and λ and setting the results to zero: We get: Sensors 2015, 15,[30827][30828][30829][30830][30831][30832][30833][30834][30835][30836][30837][30838] We putμ andλ into Equation (7). The logarithm of the profile likelihood function is expressed as: The MLE of β 0 can be obtained by differentiating Equation (9) with respect to β 0 and setting the results to zero: So far, we managed to use T 1,i , T 2,i , T 3,i , T 4,i ( N i"1 to estimate the clock offset β 0 between the two nanomachines. The proof of the convergence of the estimated offset is shown in Equation (11). The estimated offset converges in probability to the true offset value between the transmitter and the receiver: lim nÑ8`β 0´β0" lim nÑ8˜ř N i"1`T3,i`T2,i´T4,i´T1,i2 The estimator is unbiased as proved below. Equations (5) and (6) are used in the derivation:

Cramer Rao Lower Bound
In this section, the Cramer Rao Lower Bound (CRLB) is derived. The CRLB is a bound of the variance of any unbiased estimator. It is critical to evaluate how much performance loss is incurred by the proposed estimator of the clock offset based on the observations. For the molecular transmission Sensors 2015, 15,[30827][30828][30829][30830][30831][30832][30833][30834][30835][30836][30837][30838] from node A to node B in the ith round of message exchange, the likelihood function of the unknown β 0 is written as: where µ and λ are unknown, and T 1, i and T 2, i are random variables. The score is defined to indicate how sensitively a likelihood function depends on its parameter. The score for β 0 can be expressed as: The variance of the score is known as the Fisher information. It can be expressed as: (   15) where X i " T 2, i´T1, i´β0 . The closed-form expression of the expected value is difficult to compute. FIM X i pβ 0 q can be approximated as: For the message transmission from node B to node A, the Fisher information FIM Y i pβ 0 q is calculated. The result is that FIM Y i pβ 0 q is equal to FIM X i pβ 0 q. For the total N round bi-directional message exchanges, the fisher information is expressed as: The CRLB is the inverse of the fisher information. The variance of the estimator will be greater than or equal to the CRLB as:

Simulation Results and Discussions
To evaluate the proposed maximum-likelihood estimator for the clock offset, extensive simulation experiments using MATLAB have been conducted. Suppose that there are two nanomachines, a transmitter nanomachine and a receiver nanomachine. They are separated by a distance of 1-10 µm, which is a typical distance used in [37,38]. The medium velocity is set to 1-5 µm/µs [21,39]. Using the equation µ " d v mentioned in Section 2, we choose the average propagation delay of the information molecules from 1 to 5 µs which is within the feasible range. The diffusion coefficient D is set to 1-10 µm 2 /µs, which is used in [36,40]. According to the equation mentioned in Section 2, we choose the shape parameter λ " 1, which is within the feasible range. Assume that at each integral time instant, the transmitter nanomachine sends a timing-request command to the receiver nanomachine. The propagation delays of the message transmissions between the transmitter nanomachine and the receiver nanomachine are generated based on the inverse Gaussian distribution with mean µ and the shape parameter λ. Since the simulation aims to evaluate the performance of the proposed maximum likelihood estimator for the clock offset based on the system model in Equations (2) and (3), we generate the propagation delay directly based on the inverse Gaussian distribution. The initial clock offset in the simulation is set to 100 µs to 300 µs. The simulations are performed 100 times. The mean square error (MSE) of the clock offsets is calculated. Figure 3 shows the relationship between the MSE of the estimated offset and the number of rounds for different pre-defined β 0 values. It is clear that as the number of synchronization rounds increases, the MSE of the estimated clock offset decreases quickly. This result illustrates the effectiveness of our proposed estimation algorithm. It should also be noted that after 50 rounds of the synchronization, the MSE still does not reach zero. In the figure, the MSEs of the estimated clock offset for three different initial clock offset values, 100, 200 and 300 µs, are presented. The result implies that the initial clock offset does not affect the clock synchronization accuracy by the proposed clock synchronization algorithm. time instant, the transmitter nanomachine sends a timing-request command to the receiver nanomachine. The propagation delays of the message transmissions between the transmitter nanomachine and the receiver nanomachine are generated based on the inverse Gaussian distribution with mean μ and the shape parameter λ. Since the simulation aims to evaluate the performance of the proposed maximum likelihood estimator for the clock offset based on the system model in Equations (2) and (3), we generate the propagation delay directly based on the inverse Gaussian distribution. The initial clock offset in the simulation is set to 100 μs to 300 μs. The simulations are performed 100 times. The mean square error (MSE) of the clock offsets is calculated. Figure 3 shows the relationship between the MSE of the estimated offset and the number of rounds for different pre-defined 0 values. It is clear that as the number of synchronization rounds increases, the MSE of the estimated clock offset decreases quickly. This result illustrates the effectiveness of our proposed estimation algorithm. It should also be noted that after 50 rounds of the synchronization, the MSE still does not reach zero. In the figure, the MSEs of the estimated clock offset for three different initial clock offset values, 100, 200 and 300 μs, are presented. The result implies that the initial clock offset does not affect the clock synchronization accuracy by the proposed clock synchronization algorithm.   Figure 4 shows the relationship between the MSE of the estimated clock offset and the number of clock synchronization rounds with different means of the propagation delay µ. Like Figure 3, the MSE of the estimated clock offset decreases as the number of synchronization rounds increases for all the different propagation delays. This also proves the effectiveness of our proposed algorithm. It is obvious that the larger the mean propagation delay is, the larger the MSE of the estimated clock offset is. The reason for this fact is that if the mean of the propagation delay is larger, then there will be a longer period for the information molecules to move randomly in the environment. The delay will be more diverse and decrease the clock synchronization accuracy. Generally speaking, a larger mean propagation delay is due to a longer distance between the nanomachines. Therefore, a densely deployed nanonetwork could achieve a better clock synchronization performance than a loosely deployed nanonetwork in terms of accuracy. more diverse and decrease the clock synchronization accuracy. Generally speaking, a larger mean propagation delay is due to a longer distance between the nanomachines. Therefore, a densely deployed nanonetwork could achieve a better clock synchronization performance than a loosely deployed nanonetwork in terms of accuracy. It is clear that the CRLB becomes smaller as the number of message exchange rounds increases. This can be explained by Equation (16). Due to the additivity property of the Fisher information [41], as long as the experiments are independent, every observation contains the information about the unknown parameter. Therefore, as N increases, more information about β 0 is obtained. Then the CRLB becomes smaller. From Figure 4, it can also be seen that the MSE of the estimated clock offset is bigger than the CRLB. This implies that the proposed estimator is not the minimum variance unbiased (MVU) estimator. The work can be continued to find the MVU estimator of the clock offset. Figure 5 compares the MSEs of the estimated clock offsets by using the proposed MLE algorithm, the direct synchronization, and the estimator proposed in [42]. In the direct synchronization algorithm, the receiver ignores the propagation delay, and considers that the difference between its recorded receiving time stamp and the received transmitting time stamp from the transmitter is the clock offset. The algorithm in [42] is similar to the algorithm proposed in this paper but considers both the clock offset and the clock skew estimation. For fair comparison, all the three algorithms use the same input data in the simulations. From the figure, it can be seen that as the number of rounds increases, the estimated offsets become smaller, which manifests the effectiveness of all three algorithms. The MSEs It is clear that the CRLB becomes smaller as the number of message exchange rounds increases. This can be explained by Equation (16). Due to the additivity property of the Fisher information [41], as long as the experiments are independent, every observation contains the information about the unknown parameter. Therefore, as N increases, more information about β 0 is obtained. Then the CRLB becomes smaller. From Figure 4, it can also be seen that the MSE of the estimated clock offset is bigger than the CRLB. This implies that the proposed estimator is not the minimum variance unbiased (MVU) estimator. The work can be continued to find the MVU estimator of the clock offset. Figure 5 compares the MSEs of the estimated clock offsets by using the proposed MLE algorithm, the direct synchronization, and the estimator proposed in [42]. In the direct synchronization algorithm, the receiver ignores the propagation delay, and considers that the difference between its recorded receiving time stamp and the received transmitting time stamp from the transmitter is the clock offset. The algorithm in [42] is similar to the algorithm proposed in this paper but considers both the clock offset and the clock skew estimation. For fair comparison, all the three algorithms use the same input data in the simulations. From the figure, it can be seen that as the number of rounds increases, the estimated offsets become smaller, which manifests the effectiveness of all three algorithms. The MSEs of the estimated clock offset by the proposed MLE scheme is smaller than that by the other two schemes. This proves that the proposed MLE can achieve a better accuracy performance. The reason that the direct synchronization has the biggest MSE is because the propagation delay is neglected in the calculation. The model is simplified compared to the other two algorithms, while in a real scenario, the molecular communication incurs a propagation delay, which affects the estimation significantly. The reason that the proposed algorithm outperforms the estimator described in [42] is because the estimation in [42] is a joint estimation for the clock offset and the clock skew. The result is obtained when both of the two estimators achieve their maximum likelihood. For the estimation of the clock offset in this paper, the proposed estimator is derived only for the clock offset itself achieving the maximum likelihood. Therefore, the proposed solution gives better performance. timation in [42] is a joint estimation for the clock offset and the clock skew. The result is obtaine en both of the two estimators achieve their maximum likelihood. For the estimation of the cloc fset in this paper, the proposed estimator is derived only for the clock offset itself achieving th ximum likelihood. Therefore, the proposed solution gives better performance.

Figure 5.
Comparison among the proposed algorithm, direct synchronization, and estimator proposed in [42]; β 0 = 100 μs, μ = 5 μs and λ = 1. Table 1 summarizes the computational complexity of the three estimators. The numbers of th erations for the clock synchronization by the nanomachine or sensor node are shown in the tabl e values are calculated based on the expressions of the different estimators. From the figure it ar that the direct synchronization method uses fewer operations than the other two estimators. B e performance in terms of the accuracy is the worst, which is obtained in Figure 5. The complexity r proposed algorithm and the estimator proposed in [42] are both ( ).

Conclusions
In this paper, we have addressed the issue of clock synchronization among nanomachines onanosensor networks. A two-way packet-based synchronization method has been proposed. Th Figure 5. Comparison among the proposed algorithm, direct synchronization, and estimator proposed in [42]; β 0 = 100 µs, µ = 5 µs and λ = 1. Table 1 summarizes the computational complexity of the three estimators. The numbers of the operations for the clock synchronization by the nanomachine or sensor node are shown in the table. The values are calculated based on the expressions of the different estimators. From the figure it is clear that the direct synchronization method uses fewer operations than the other two estimators. But the performance in terms of the accuracy is the worst, which is obtained in Figure 5. The complexity of our proposed algorithm and the estimator proposed in [42] are both O pNq.

Conclusions
In this paper, we have addressed the issue of clock synchronization among nanomachines in bionanosensor networks. A two-way packet-based synchronization method has been proposed. The maximum-likelihood estimator for the clock offset between nanomachines has been derived. The simulation results have demonstrated the effectiveness of the proposed scheme. The convergence and the bias of the estimator have also been verified. The work presented in this paper is a very important step towards achieving distributed cooperation among nanomachines in bionanosensor networks. Future work will focus on solutions to achieve global clock synchronization in nanosensor networks.