Membrane Fusion-Based Transmitter Design for Static and Diffusive Mobile Molecular Communication Systems

We assume that the spherical TX releases molecules from its outer membrane after fusion between the membrane and vesicles. Each vesicle stores and transports $\eta$ type $\sigma$ molecules. We consider that the TX is filled with the fluid medium that has uniform temperature and viscosity. Generally, vesicles can diffuse in the fluid medium based on the experiment in [19] or can be actively transported along microtubules by motor proteins [20]. In this paper, we focus on the diffusion of vesicles. We assume that a vesicle cluster [21] exists within the TX. Vesicle clustering is a biological phenomenon observed in the presynaptic neuron, where vesicles can form clusters to sustain a high release probability of neurotransmitters. Vesicles can be released from the cluster via dephosphorylation [22]. Since we assume that the size of a vesicle [23] is relatively small compared to the size of the TX, we mathematically model the vesicle release process as an instantaneous release from a point within the TX. Once vesicles are released, they diffuse randomly with a constant diffusion coefficient $D_{\mathrm{v}}$. According to [9], the natural fusion of a vesicle and the cell membrane can be considered as two steps: 1) v-SNARE proteins ($\mathrm{S}_{\mathrm{v}}$) on the vesicle membrane bind to t-SNARE proteins ($\mathrm{S}_{\mathrm{t}}$) on the cell membrane to generate trans-SNARE complexes ($\mathrm{S}_{\mathrm{c}}$), and 2) $\mathrm{S}_{\mathrm{c}}$ catalyzes the fusion of vesicular and cell membranes. For tractability, we apply the following main assumptions to the TX model:

• A1)

  Vesicles are released from the TX’s center. This assumption ensures that molecules are uniformly released from the TX membrane. Considering vesicles released from any point within the TX is an interesting topic for future work.

• A2)

  The interaction between $\mathrm{S}_{\mathrm{v}}$ and $\mathrm{S}_{\mathrm{t}}$ is modeled as the irreversible reaction $\mathrm{S}_{\mathrm{v}}+\mathrm{S}_{\mathrm{t}}\rightarrow\mathrm{S}_{\mathrm{c}}$, which indicates that we only focus on the forward reaction between $\mathrm{S}_{\mathrm{v}}$ and $\mathrm{S}_{\mathrm{t}}$. The backward reaction for the disassembly of $\mathrm{S}_{\mathrm{c}}$ occurs after MF, wherein the free $\mathrm{S}_{\mathrm{v}}$ and $\mathrm{S}_{\mathrm{t}}$ are used for the subsequent rounds of transport [24]. This is beyond the scope of our system model. After $\mathrm{S}_{\mathrm{c}}$ is generated, it catalyzes MF, which is $\mathrm{S}_{\mathrm{c}}\rightarrow\mathrm{MF}$ [25]. We simplify this reaction network by removing the intermediate product $\mathrm{S}_{\mathrm{c}}$ and combining these two reactions [26], which is given by

  $\displaystyle\mathrm{S}_{\mathrm{v}}+\mathrm{S}_{\mathrm{t}}\stackrel{{\scriptstyle k_{\mathrm{f}}}}{{\longrightarrow}}\mathrm{MF},$

  (1)

  where $k_{\mathrm{f}}$ is the forward reaction rate in $\mu\mathrm{m}/\mathrm{s}$.

• A3)

  The membrane is fully covered by an infinite number of $\mathrm{S}_{\mathrm{t}}$ and the occupancy of $\mathrm{S}_{\mathrm{t}}$ by $\mathrm{S}_{\mathrm{v}}$ is ignored. Assuming perfect receptor coverage and ignoring occupancy are for the sake of tractability and have been widely adopted in previous studies, e.g., [27, 28].

• A4)

  Once molecules are released, the spherical TX is assumed to not hinder the random diffusion of molecules in the propagation environment, i.e, the TX is transparent to the diffusion of released molecules. [4] analyzed the hindrance of the TX membrane to the diffusion of molecules using simulation. The joint consideration of the hindrance of the TX membrane and the absorbing RX is cumbersome for theoretical analysis. We will investigate the impact of the TX membrane on molecule diffusion in the propagation environment in future work. In Fig. 7(b) and Table IV, we relax this assumption by treating the TX membrane as a reflecting boundary for molecules in the propagation environment and perform the corresponding particle-based simulation [29]. Results indicate that an obvious deviation is caused only when the TX and RX are very close to each other.

Based on assumptions A2 and A3, if a vesicle hits the TX membrane, it fuses to the membrane with a probability of $k_{\mathrm{f}}\sqrt{\frac{\pi\Delta t}{D_{\mathrm{v}}}}$ during a time interval $\Delta t$ [30]. We define this probability as the MF probability²²2Based on this definition, the MF-based TX can be regarded as a TX with a semi-permeable boundary as in [6]. We clarify that the major difference between our work and [6] is the method of analysis. The authors in [6] applied the transfer function approach to investigate the molecule concentration within the $\mathrm{TX}$, while we focus on the end-to-end communication channel and derive analytical expressions for the CIR in this channel.. The equation of the MF probability is accurate if $k_{\mathrm{f}}\ll\sqrt{\frac{D_{\mathrm{v}}}{\pi\Delta t}}$. The MF probability is applied in the simulation process in Section VI. After MF, the type $\sigma$ molecules stored by the vesicle are instantaneously released into the propagation environment. The location and time for the occurrence of MF are the initial location and time for molecules to start moving in the propagation environment. We note that the time scale for the MF process is ignored since it is relatively small compared to the entire transmission process [31].

We consider the following two cases for TX and RX mobility:

We assume for both static and mobile scenarios that information transmitted from the TX to the RX is encoded into a sequence of binary bits with length $W$. The sequence of binary bits is $\mathbf{b}_{1:W}=\left[b_{1},b_{2},...b_{W}\right]$ and $b_{w}$ is the $w$th bit. We denote $\phi$ as the bit interval length. We assume that each bit is transmitted at the beginning of the current bit interval with probabilities $\mathrm{Pr}(b_{w}=1)=P_{1}$ and $\mathrm{Pr}(b_{w}=0)=P_{0}=1-P_{1}$, where $\mathrm{Pr}(\cdot)$ denotes the probability. For the static scenario, we consider that $b_{1}$ is transmitted at time $t=0$. For the mobile scenario, we consider that $b_{1}$ is transmitted at time $t^{\prime}$ such that the distance between the centers of the TX and RX is a random variable (RV) when $b_{1}$ is transmitted. We adopt ON/OFF keying for the modulation, which means that $N_{\mathrm{v}}$ vesicles are released from the TX’s center at the beginning of the bit interval to transmit bit $1$ while no vesicle is released to transmit bit 0. We further assume that the TX and RX are perfectly synchronized since several practical synchronization schemes for MC have been studied and proposed in [32]. For demodulation of $b_{w}$ at the RX, we adopt a widely-investigated threshold detector that compares the number of molecules absorbed during the $w$th bit interval with a detection threshold [16]. In this paper, we consider the threshold of the detector is fixed for simplicity. More complex detectors, e.g., an adaptive threshold detector, are summarized in [33].

As MF guarantees the release of molecules from vesicles to the propagation environment, the probability of releasing molecules equals the vesicle fusion probability. Thus, we need to obtain the distribution of vesicle concentration within the TX to derive the release probability from the TX membrane. In the spherical coordinate system, we denote $C(r,t)$, $0\leq r\leq r_{\scriptscriptstyle\textnormal{T}}$, as the distribution of vesicle concentration at time $t$ and distance $r$ from the TX’s center. When an impulse of vesicles is released from the TX’s center at $t=0$, the initial condition is expressed as [34, eq. 3(c)]

where $\delta(\cdot)$ is the Dirac delta function. We then apply Fick’s second equation to describe the diffusion of vesicles within the TX. In general, Fick’s second equation is a spatially 3D partial differential equation. As we consider that vesicles are released from the TX’s center, the TX model is spherically symmetric such that the diffusion of vesicles is described as [35, eq. (2.7)]

Based on assumption A2, the boundary condition is described by the irreversible reaction given by (1), which can be characterized by the radiation boundary condition [36, eq. (3.25)]

where the left-hand side represents the diffusion flux³³3The diffusion flux $[\mathrm{molecule}\cdot\mathrm{m}^{-2}\cdot\mathrm{s}^{-1}]$ is the rate of movement of molecules across a unit area in unit time [35]. over the TX membrane. Since the flux is in the positive radial direction, the right-hand side is positive.

negative sign on the right-hand side indicates that the condition is over the inner boundary.

Based on the initial condition in (2), Fick’s second law in (3), and the boundary condition in (4), we derive the analytical expression for the molecule release probability from the TX membrane, denoted by $f_{\mathrm{r}}(t)$, in the following theorem:

We denote $F_{\mathrm{r}}(t)$ as the fraction of molecules released by time $t$ and obtain it via $F_{\mathrm{r}}(t)=\int_{0}^{t}f_{\mathrm{r}}(u)\mathrm{d}u$. We present $F_{\mathrm{r}}(t)$ in the following corollary:

As each molecule has the same probability to be released by time $t$, i.e., $F_{\mathrm{r}}(t)$, the number of molecules released by time $t$ is $N_{\mathrm{v}}\eta F_{\mathrm{r}}(t)$.

In this subsection, we derive the hitting probability when molecules are uniformly released from the TX membrane, i.e., ignoring the internal molecules’ propagation within vesicles inside the TX and the TX’s MF process. This hitting probability establishes the foundation for deriving the end-to-end molecule hitting probability at the RX surface in Section III-C. We consider that molecules are initially uniformly distributed over the TX membrane and released simultaneously at $t=0$. We denote $\Omega_{\scriptscriptstyle\textnormal{T}}$ as the membrane area and $p_{\mathrm{u}}(t)$ as the hitting probability at the RX due to the uniform release of molecules. We clarify that uniformly-distributed molecules means that the likelihood of a molecule released from any point on the TX membrane is the same. We denote this probability by $\rho$, where we have $\rho=\left(4\pi r_{\scriptscriptstyle\textnormal{T}}^{2}\right)^{-1}$. We further consider an arbitrary point $\alpha$ on the TX membrane. Based on [38, eq. (9)], the hitting probability of a molecule at the RX at time $t$ when the molecule is released from the point $\alpha$ at time $t=0$, $p_{\alpha}(t)$, is given by

where $l_{\alpha}$ is the distance between the point $\alpha$ and the center of the RX.

Given that molecules are distributed uniformly over the TX membrane, $p_{\mathrm{u}}(t)$ is obtained by taking the surface integral of $p_{\alpha}(t)$ over the spherical TX membrane. Using this method, we solve $p_{\mathrm{u}}(t)$ in the following lemma:

We denote $P_{\mathrm{u}}(t)$ as the fraction of molecules absorbed by time $t$ when the TX uniformly releases molecules and obtain it via $P_{\mathrm{u}}(t)=\int_{0}^{t}p_{\mathrm{u}}(u)\mathrm{d}u$. We present $P_{\mathrm{u}}(t)$ in the following corollary:

We denote $p_{\mathrm{e}}(t)$ as the end-to-end molecule hitting probability at the RX when an impulse of vesicles is released from the TX’s center at time $t=0$. The molecule release probability from the TX membrane at time $u$, $0\leq u\leq t$, is given by (5), which is $f_{\mathrm{r}}(u)$. For molecules released at time $u$, the molecule hitting probability at the RX at time $t$ is given by (9), which is $p_{\mathrm{u}}(t-u)$. Based on $f_{\mathrm{r}}(u)$ and $p_{\mathrm{u}}(t-u)$, $p_{\mathrm{e}}(t)$ is given by $p_{\mathrm{e}}(t)=\int_{0}^{t}p_{\mathrm{u}}(t-u)f_{\mathrm{r}}(u)\mathrm{d}u$. Applying (5) and (9) to this equation, we derive $p_{\mathrm{e}}(t)$ in the following theorem:

We denote $P_{\mathrm{e}}(t)$ as the end-to-end fraction of molecules absorbed at the RX by time $t$ and obtain it via $P_{\mathrm{e}}(t)=\int_{0}^{t}p_{\mathrm{e}}(u)\mathrm{d}u=\int_{0}^{t}P_{\mathrm{u}}(t-u)f_{\mathrm{r}}(u)\mathrm{d}u$. We present $P_{\mathrm{e}}(t)$ in the following corollary:

We further denote $P_{\mathrm{e},\infty}$ as the end-to-end asymptotic fraction of molecules absorbed as $t\rightarrow\infty$. We present $P_{\mathrm{e},\infty}$ in the following corollary:

We denote $g(l_{t^{\prime}}|l_{0})$ as the probability density function (PDF) of $l_{t^{\prime}}$, which is the distance between centers of the TX and RX at time $t^{\prime}$, given that the initial distance between the centers of the TX and RX is $l_{0}$ at $t=0$. The coordinates of the TX’s center and RX’s center at time $t^{\prime}$ are $\left[x_{\scriptscriptstyle\textnormal{T}}(t^{\prime}),y_{\scriptscriptstyle\textnormal{T}}(t^{\prime}),z_{\scriptscriptstyle\textnormal{T}}(t^{\prime})\right]$ and $\left[x_{\scriptscriptstyle\textnormal{R}}(t^{\prime}),y_{\scriptscriptstyle\textnormal{R}}(t^{\prime}),z_{\scriptscriptstyle\textnormal{R}}(t^{\prime})\right]$, respectively. As both TX and RX perform Brownian motion, the displacement of each coordinate follows a Gaussian distribution. $l_{t^{\prime}}$ is calculated as $l_{t^{\prime}}=\sqrt{(x_{\scriptscriptstyle\textnormal{T}}(t^{\prime})-x_{\scriptscriptstyle\textnormal{R}}(t^{\prime}))^{2}+(y_{\scriptscriptstyle\textnormal{T}}(t^{\prime})-y_{\scriptscriptstyle\textnormal{R}}(t^{\prime}))^{2}+(z_{\scriptscriptstyle\textnormal{T}}(t^{\prime})-z_{\scriptscriptstyle\textnormal{R}}(t^{\prime}))^{2}}$, where $x_{\scriptscriptstyle\textnormal{T}}(t^{\prime})-x_{\scriptscriptstyle\textnormal{R}}(t^{\prime})$, $y_{\scriptscriptstyle\textnormal{T}}(t^{\prime})-y_{\scriptscriptstyle\textnormal{R}}(t^{\prime})$, and $z_{\scriptscriptstyle\textnormal{T}}(t^{\prime})-z_{\scriptscriptstyle\textnormal{R}}(t^{\prime})$ follow the Gaussian distribution. Thus, $\frac{l_{t^{\prime}}}{\sqrt{2D_{1}t^{\prime}}}$ follows a noncentral chi distribution with three degrees of freedom [39], and $g(l_{t^{\prime}}|l_{0})$ is given by [16, eq. (8)]

where $D_{1}=D_{\scriptscriptstyle\textnormal{T}}+D_{\scriptscriptstyle\textnormal{R}}$ and $I_{\frac{1}{2}}(\cdot)$ is the modified Bessel function of the first kind [36]. We note that $D_{\scriptscriptstyle\textnormal{T}}$, $D_{\scriptscriptstyle\textnormal{R}}$, and $t^{\prime}$ should not be very large. Otherwise, $l_{t^{\prime}}$ may become less than $r_{\scriptscriptstyle\textnormal{T}}+r_{\scriptscriptstyle\textnormal{R}}$ such that the TX and RX overlap with each other.

We denote $p_{\mathrm{u,m}}(\tau,t^{\prime})$ as the expected hitting probability at the mobile RX at time $t^{\prime}+\tau$ when the mobile TX releases molecules uniformly over the TX membrane at time $t^{\prime}$. To describe the relative motion between released molecules and the RX, we define an effective diffusion coefficient $D_{2}$ as $D_{2}=D_{\scriptscriptstyle\textnormal{R}}+D_{\sigma}$ [40]. By replacing $l$, $t$, and $D_{\sigma}$ in (9) with $l_{t^{\prime}}$, $\tau$, and $D_{2}$, respectively, we obtain the hitting probability at time $t^{\prime}+\tau$, i.e., $p_{\mathrm{u}}(\tau)\big{|}_{l=l_{t^{\prime}},D_{\sigma}=D_{2}}$, for a given $l_{t^{\prime}}$. We note that $l_{t^{\prime}}$ is a RV with the PDF given by (15). Hence, we obtain $p_{\mathrm{u,m}}(\tau,t^{\prime})$ as $p_{\mathrm{u,m}}(\tau,t^{\prime})=\int_{0}^{\infty}g(l_{t^{\prime}}|l_{0})p_{\mathrm{u}}(\tau)\big{|}_{l=l_{t^{\prime}},D_{\sigma}=D_{2}}\mathrm{d}l_{t^{\prime}}$ and present it in the following lemma:

We denote $p_{\mathrm{e,m}}(\tau,t^{\prime})$ as the expected end-to-end molecule hitting probability at the mobile RX at time $t^{\prime}+\tau$ when an impulse of vesicles is released from the TX’s center at time $t^{\prime}$. The molecule release probability from the TX at time $t^{\prime}+u$, $0\leq u\leq\tau$, is given by (5), which is $f_{\mathrm{r}}(u)$. For molecules released at time $t^{\prime}+u$, the expected hitting probability at time $t^{\prime}+\tau$ at the RX is given by (2), which is $p_{\mathrm{u,m}}(\tau-u,t^{\prime}+u)$. Based on that, $p_{\mathrm{e,m}}(\tau,t^{\prime})$ is given by $p_{\mathrm{e,m}}(\tau,t^{\prime})=\int_{0}^{\tau}f_{\mathrm{r}}(u)p_{\mathrm{u,m}}(\tau-u,t^{\prime}+u)\mathrm{d}u$. Applying (5) and (2) to this equation, we derive $p_{\mathrm{e,m}}(\tau,t^{\prime})$ in the following theorem:

Due to the delay of molecule transport, the RX may receive molecules released from the current and all previous bit intervals. According to (3), the fraction of molecules absorbed during the $w$th bit interval when vesicles for storing those molecules released at the beginning of the $i$th bit interval, $0\leq i\leq w$, is expressed as $P_{\mathrm{e}}((w-i+1)\phi)-P_{\mathrm{e}}((w-i)\phi)$, where $\phi$ is the bit interval length. We denote $N_{w,i}$ as the number of molecules absorbed during the $w$th bit interval for transmitting $b_{i}$. As the diffusion of vesicles and released molecules are independent and molecules have the same probability to be absorbed, $N_{w,i}$ can be approximated as a Poisson RV when the number of emitted vesicles is large and the success probability $P_{\mathrm{e}}(t)$ is small [41], i.e., $N_{w,i}\sim b_{i}\mathrm{Poiss}\left(N_{\mathrm{v}}\eta\left(P_{\mathrm{e}}((w-i+1)\phi)-P_{\mathrm{e}}((w-i)\phi)\right)\right)$. Here, $\mathrm{Poiss}(\cdot)$ refers to a Poisson distribution. This is because sufficient diversity in the vesicle fusion points over the TX membrane is required to match the derivation of the end-to-end channel that includes integration over the entire TX membrane. As the sum of independent Poisson RVs is also a Poisson RV, we model the total number of molecules absorbed during the $w$th bit interval, denoted by $N_{w}$, as

where $\psi=N_{\mathrm{v}}\eta\sum_{i=1}^{w}b_{i}\left(P_{\mathrm{e}}((w-i+1)\phi)-P_{\mathrm{e}}((w-i)\phi)\right)$. Based on (19), the cumulative distribution function (CDF) of the Poisson RV $N_{w}$ is written as

In Table I, we calculate the coefficient of determination [42], denoted by $R^{2}$, for the CDF of $N_{w}$, where we vary $N_{\mathrm{v}}$ and $\eta$ and keep the total number of molecules emitted as 1000. We set $w=1$, $D_{\mathrm{v}}=18\;\mu\mathrm{m}^{2}/\mathrm{s}$, $k_{\mathrm{f}}=30\;\mu\mathrm{m}/\mathrm{s}$, and $b_{w}=1$. For other parameter values, please see Table III. For the simulation details, please see Section VI. The $R^{2}$ is used to measure the preciseness between simulation and theoretical analysis. Specifically, $R^{2}$ is given by $R^{2}=1-SS_{\mathrm{res}}/SS_{\mathrm{tot}}$, where $SS_{\mathrm{res}}$ is the sum of squared differences between the theoretical analysis and simulation, and $SS_{\mathrm{tot}}$ is the sum of squared differences between the simulation and its mean. Since $SS_{\mathrm{res}}/SS_{\mathrm{tot}}$ represents the fraction of variance unexplained, the closer the value of $R^{2}$ is to 1, the more accurate the theoretical analysis. In this table, we observe that the $R^{2}$ increases with an increase in the number of vesicles from 10 to 50, which indicates that the accuracy of (20) increases when the number of vesicles emitted increases. This is due to the fact that the Poisson approximation is impacted by two sets of trials, i.e., the number of emitted vesicles and molecules, and the approximation is more sensitive to the number of emitted vesicles since it determines the initialization of the emission of molecules. When $N_{\mathrm{v}}$ is larger than 50, we observe a small variation of $R^{2}$ for $N_{\mathrm{v}}$ from 50 to 1000. We further increase the number of realizations applied in Table I by ten times for the simulation when $N_{\mathrm{v}}=10$ and $\eta=100$, and obtain that $R^{2}=0.9771$. This indicates that the number of realizations applied for the simulation in Table I is sufficient.

We adopt a simple threshold detector, where $N_{w}$ is compared with a detection threshold to demodulate $b_{w}$. We present the threshold detector as

where $\hat{b}_{w}$ is the demodulated bit of $b_{w}$ at the RX and $\xi$ is the detection threshold. Based on (20) and (23), the error probability of $\hat{b}_{w}$ for $b_{w}=1$ and $b_{w}=0$ are

and

respectively. The BER of the $w$th bit, denoted by $Q_{\textnormal{s}}[w|\mathbf{b}_{1:w-1}]$, is given by

We denote $Q_{\textnormal{s}}$ as the average BER over all realizations of $\mathbf{b}_{1:w-1}$ and all bits from $1$ to $W$ for the static scenario. We present $Q_{\textnormal{s}}$ as

where $\Psi_{w}$ stands for all realizations of $\mathbf{b}_{1:w-1}$.