Design and Optimizing of On-Chip Kinesin Substrates for Molecular Communication
11 Design and Optimizing of On-Chip KinesinSubstrates for Molecular Communication
Nariman Farsad,
Student Member, IEEE,
Andrew W. Eckford,
Member, IEEE, and Satoshi Hiyama,
Member, IEEE,
Abstract
Lab-on-chip devices and point-of-care diagnostic chip devices are composed of many differentcomponents such as nanosensors that must be able to communicate with other components withinthe device. Molecular communication is a promising solution for on-chip communication. In particular,kinesin driven microtubule (MT) motility is an effective means of transferring information particles fromone component to another. However, finding an optimal shape for these channels can be challenging.In this paper we derive a mathematical optimization model that can be used to find the optimal channelshape and dimensions for any transmission period. We derive three specific models for the rectangularchannels, regular polygonal channels, and regular polygonal ring channels. We show that the optimalchannel shapes are the square-shaped channel for the rectangular channel, and circular-shaped channelfor the other classes of shapes. Finally, we show that among all 2 dimensional shapes the optimal designchoice that maximizes information rate is the circular-shaped channel.
Index Terms
Microfluidic Channels, Molecular Communication, Optimal Channel Design, Channel Capacity,Kinesin Substrate, Microtubule Motility, Active Transport.
This work as been presented in part at 2012 IEEE NANO conference.Nariman Farsad and Andrew W. Eckford are with the Department of Electrical Engineering and Computer Science, YorkUniversity, 4700 Keele Street, Toronto, Ontario, Canada M3J 1P3. Emails: [email protected], [email protected] Hiyama is with Research Laboratories, NTT DOCOMO Inc., Yokosuka, Kanagawa, Japan. Email: [email protected] a r X i v : . [ c s . I T ] M a y I. I
NTRODUCTION
With the advancements in the field of nanotechnology, applications such as lab-on-chip devices[1], and point-of-care diagnostic chips [2] are becoming a reality. In most of these applicationscommunication between different components in each device, or between a biological entity anda component in the device, such as a nanosensor, is required. Therefore, one of the obstaclesthat must be overcome before many of these applications can be fully realized, is devising acommunication system for small scales [3], [4].Inspired by nature, one of the most promising solutions is molecular communication [5], [6],where molecules released by a transmitter device propagate in a fluidic environment and transferinformation to a receiver device. Information can be conveyed by encoding messages into therelease timing [7], number [8], [9], concentration [10], or identities of particles [11]. Differenttypes of propagation are possible for transporting information particles in on-chip molecularcommunication such as: diffusion [12]–[14], diffusion with flow [15], active transport usingmolecular motors and cytoskeletal filaments [16], bacterial assisted propagation [11], [17], [18],and kinesin molecular motors moving over immobilized microtubule (MT) tracks [19], [20].Molecular communication is biocompatible, and can be very energy efficient. Because of theseproperties, molecular communication has attracted a lot of attention in recent years [21]. It hasmany potential applications in biomedical engineering from lab-on-chip devices [22] and point-of-care diagnostic devices [23] to targeted drug delivery [24]. However, most previous workshave considered diffusion based molecular communication, which can be very slow for on-chipapplications. Moreover, in [9] it was shown that active transport using MT filaments glidingover kinesin covered substrate can generate higher information rates compared to flow basedpropagation over larger separations distances. Other types of propagations are either difficult toimplement or have not been implemented yet. Moreover, in [25]–[27] it is shown that electricalcurrents can be used to control the speed and direction of the MTs. This makes kinesin-driven MTmotility a suitable form of propagation for on-chip molecular communication and transport, withpotential applications in lab-on-chip devices, point-of-care diagnostic devices, organ-on-a-chipdevices, and microfluid devices.For the kinesin-microtubule based molecular communication to be fully realized, guide postsfor design and implementation are required. For example, in microfluidics it is well known that different channel designs could be used to implement different functionalities in droplet microflu-idics [28]. However, finding the optimal design for these channels can be difficult because of thesystem complexities such as the random motion of the MTs. Moreover, experimental trial-and-error-based design is very laborious and time-consuming. Therefore, simulation environmentsbased on experimentally validated models such as [29], [30] can be a very effective tool forproviding design guide posts. For example in our previous works [8], [31], optimization ofthe shape of the transmission zone (i.e. location placement of transmitters) is considered usingcomputer simulations, and it was shown that optimal transmission zone is along of the walls ofthe channel. This is because usually MTs follow the channel walls [32]. Moreover, in [33], [34]we have shown that the shape of the channel can have a huge effect on information transmissionrate, and we provide basic guidepost for designing the shape of the channel.In this work, we extend the preliminary results presented in our previous works and derivean optimization model that could be used to find the optimal channel shape and dimensions.Our results can be used in two ways. First, in on-chip transport applications, our results couldbe used to maximize the transport rate. On-chip transport is a very important process in lab-on-chip devices. For example, in [35], kinesin-microtubule transport scheme is used for high-throughput molecular transport and assembly. Our results can also be used in on-chip molecularcommunication applications to maximize the channel capacity of the system. We show that usingthe optimized channel shape can increase the capacity by more than 1.5 bits per channel use andin some cases by more than 2.5 bits per channel use. This can be a significant gain in sequentialdata transmission.To setup our optimization formula, first we derive a complete model relating the number ofMT trips during a single transmission period and the shape of the channel. This is our first majorcontribution. We then use this model to find the optimal shape and dimensions for rectangularchannels, regular polygonal channels, and regular polygonal ring-shaped channels. This is oursecond major contribution. Finally, we show that for all possible two dimensional channel shapes,the optimal channel is the circular-shaped channel. Moreover, for a particular transmission ratewe show that we could use our optimization model to find the optimal radius for this circle.This is our third major contribution. We verify all our results using computer simulations.The rest of this paper is organized as follows. In Section II an overview of the kinesin-microtubule based molecular communication system is presented. A general information trans-
Fig. 1: Depiction of the communication environment.mission model is presented in Section III. In Section IV we derive the channel shape optimizationmodel, and in Section V we use the derived model to find the optimal channel design strategies.The results are then compared with computer simulations in Section VI. Concluding remarksare given in Section VII.II. K
INESIN -M ICROTUBULE BASED M OLECULAR C OMMUNICATION C HANNELS
For our model we assume that an area in the kinesin covered channel is designated as the transmission zone , where the transmitter and the information particles released by the transmitterreside. Similarly, an area in the channel is designated as the receiver zone , where the receiver andits receptors are placed. Generally, we assume the transmission zone and the receiver zone arelocated at the opposite ends of the channel. Moreover, the information particles released by thetransmitter remain stationary at the transmission zone until they are picked up by MT filamentsgliding over kinesin covered substrate. The MT filaments then carry these information particlesto the receiver zone, where they are unloaded and delivered to the receiver. The information canbe encoded in the type, number, release timing, and concentration of the particles released anddelivered. Moreover, the information particles can be encapsulated inside liposomes to protectthem from denaturation (e.g., molecular deformation and cleavage caused by enzymatic attacksor changes in pH of the outer aqueous phase) in the propagation environment.For the loading and unloading information particles, we assume deoxyribonucleic acid (DNA)hybridization bonds are employed as shown in [36]. In this scheme, information particles, MT l w rn = 6 r i r o n = 8Rectangular Regular Polygon Polygon Ring Fig. 2: Different shape classes and their parameters.filaments, transmission zone, and the receiver zones are covered with single stranded DNAs(ssDNA). When information particles are released by the transmitter, they anchor themselvesto the transmission zone through DNA hybridization bonds until an MT filament passes at thevicinity of the anchored particle. At this point the particle is loaded on to the passing MT againthrough DNA hybridization. When a loaded MT enters this receiver zone, the loaded particlesare unloaded using yet another DNA hybridization bond. This process is summarized in Figure1 and the reader is referred to [36] for detailed explanation and laboratory demonstration of thistechnique.
A. Shape of the Channel
In this paper we assume the channel can be formed to have a wide variety of different shapes,and we find the optimal shape that maximizes information transmission rate among all possibleshapes. In particular, we assume that channel shape can belong to any of the following threeclasses of shapes: rectangular, regular polygon, and regular polygon ring. Regular polygons, areequiangular (all angles are equal in measure) and equilateral (all sides have the same length),and they include a large class of geometric shapes, ranging from equilateral triangles, squares,pentagons, and hexagons, all the way up to a circle as number of sides approaches infinity.Figure 2 shows an example from each of the three different shape classes. In practice, thesechannels, as well as channel having any other shape can be created using similar proceduresused in [16], [30], [36].Depending on the class of the channel shape, different parameters can be used to further adjust the shape within the shape class. For rectangular channels two parameters can be used toadjust the dimensions and the shape: width w and length l . When both the width and the lengthare equal in value the channel would be a square channel. Similarly, for regular polygons twoparameters can be used to further adjust the shape: the number of sides n , and the radius of thecircumscribed circle r (this radius can also be defined as the distance from the center point of theregular polygon to one of the edges of the polygon). For ring-shaped channels three parametersare used to define the shape: the number of sides n , the radius of the inner circumscribed circle r i , and the radius of the outer circumscribed circle r o . Figure 2 summarizes all these parametersfor different classes of shapes.Regardless of the shape of the channel, without loss of generality, we assume the transmissionzone is always on the left side of the channel and the receiver zone is along the right side of thechannel. In [8], it was shown that the optimal transmission zone is along the walls of the channelsince MTs mostly glide very close to channel walls. This positioning increases the chance ofan MT picking up an information particle during its trips. Therefore, in this work we alwaysassume the transmission zone is along the walls of the channel regardless of the channel shape.Moreover, the minimum separation distance between the start of the transmission zone and thereceiver zone are always the same in our numerical evaluations of Section VI regardless of theshape of the channel. B. Simulating the Channel
In [16], [36] the authors show that creating kinesin-microtubule based channels described inthe previous section is possible in wet labs. However, studying the resulting molecular communi-cation channel, and solving design and optimization problems using laboratory experimentationis extremely difficult because these experiment are time consuming, laborious, and expensive.To overcome this issue, computer simulations have been used in previous works [8], [9], [30].Similarly, in this work we use Monte-Carlo techniques proposed in [29] to simulate the motionof the MTs over a kinesin covered substrate. The equations for the motion of the MT, whichwere developed and verified experimentally in [29], are given below.The motion of the MT is largely regular, although the effects of Brownian motion causerandom fluctuations. Since the MTs move only in the x – y directions, and do not move in the z direction (along the height of the channel), we consider a two-dimensional simulation of MTs for ∆ t time intervals. Given some initial position ( x , y ) at time t = 0 , for any integer k > ,the motion of the MT is given by the sequence of coordinates ( x i , y i ) for i = 1 , , . . . , k . Eachcoordinate ( x i , y i ) represents the position of the MT’s head at the end of the time t = i ∆ t ,where x i = x i − + ∆ r cos θ i , (1) y i = y i − + ∆ r sin θ i . (2)In this case, the step size ∆ r at each step is an independent and identically distributed (i.i.d.)Gaussian random variable with mean and variance E [∆ r ] = v avg ∆ t, (3) Var[∆ r ] = 2 D ∆ t, (4)where v avg is the average velocity of the MT, and D is the MT’s diffusion coefficient. The angle θ i is no longer independent from step to step: instead, for some step-to-step angular change ∆ θ ,we have that θ i = ∆ θ + θ i − . (5)Now, for each step, ∆ θ is an i.i.d. Gaussian-distributed random variable with mean and variance E [∆ θ ] = 0 , (6) Var[∆ θ ] = v avg ∆ tL p , (7)where L p is the persistence length of the MT’s trajectory. Following [29], in case of a collisionwith a boundary, we assume that the MT sets θ i so as to follow the boundary . Finally, oursimulator can model any polygon-shaped channel.We also use the grid loading mechanism proposed in [8], to simulate the information particles’loading and unloading. For loading an information particle, the MT filament must drive close tothe anchored particle. Therefore, we divide the transmission zone into a square grid, where thelength of each square in the grid is the same as the diameter of the particles. We then distributeparticles randomly and uniformly among the squares in the grid. In general, we assume thatthe MTs can load multiple particles, which we know to be possible based on lab experiments;thus, if an MT enters a square which is occupied by a particle, and it has an empty loading Grid TransmissionZone Microtubule StartsReceiver ZoneMicrotubuleLoads MicrotubuleUnloads(Delivers)
Fig. 3: An example of the simulation environment.slot available, we assume the MT loads that particle. For unloading, we assume all the loadedparticles are unloaded as soon as an MT enters the receiver zone.Figure 3 shows a sample simulation. In this figure, the channel is a regular polygon ring withparameters n = 8 , r i = 20 µ m, and r o = 25 µ m. Other parameters used for simulation are:simulation time steps of ∆ T = 0 . seconds, MT diffusion coefficient D = 2 . · − µ m /s,average speed of the MT v avg = 0 . µ m/s, and persistence length of the MT trajectory L p = 111 µ m. We also assume the size of the information particles is 1 µ m, the average length of the MTsis 10 µ m, and each MT can load up to 5 information particles in one trip from the transmissionzone to the receiver zone. These parameters are all selected based on experimental observationsof DNA covered MTs moving over a kinesin covered substrate, and the reader is referred to [9]for more details about the simulation environment. In the rest of the paper we will always usethese parameters unless it is specified otherwise.In Figure 3 the MT initially starts on top of the channel. It will then move towards thetransmission zone (blue/cyan grid) along the left wall of the channel. The MT loads the first 5information particles on its way and continues to the receiver zone where it unloads the loaded information particles. Notice that although the MT enters grid squares that contain informationparticles, after 5 loaded particles, it does not load anything else because the maximum loadcapacity is reached. III. I NFORMATION T RANSMISSION M ODEL
In molecular communication, messages can be encoded into information carrying particlesusing many different schemes. Information can be encoded into the number, concentration, ortype of the particles released. One key observation in molecular communication is that, regardlessof the encoding scheme, information is always transmitted through mass transfer (i.e. transferof particles). Therefore, in this section we present a general information transmission model formass transfer, which is independent of the information encoding scheme.Let X = { , , , · · · , x max } be the set of possible particles that could be released by thetransmitter, where x max is the maximum number of particles the transmitter can release. Let X ∈ X be the number of information particles released into the medium by the transmitter.We define time per channel use (TPCU) as a predefined amount of time T representing theduration of a single transmission session. Because of the random motion of the MT filaments,given this time limit, information might not be perfectly conveyed to the receiver. For example,it is possible that some of the particles will not arrive at the receiver after T has elapsed andtherefore there is some information loss. This effect is similar in nature to the noise introducedby conventional electronic or electromagnetic channels.Let Y ( T ) ∈ X represent the number of information particles that arrive at the destinationafter TPCU duration T . Because the receiver does not have any prior knowledge of the numberof particles released by the transmitter, at the receiver X is a discrete random variable givenby probability mass function (PMF) f X ( x ) . Similarly, given X particles were released by thetransmitter, and the TPCU duration T , Y ( T ) is a discrete random variable given by conditionalPMF f Y ( T ) | X ( y ( T ) | x ) . This conditional PMF is very important because it characterizes thechannel completely. It can be used to calculate important channel parameters such as channelcapacity [37], the maximum rate at which any communication system can reliably transmitinformation over a noisy channel.Assuming that the transmitter and receiver are perfect (i.e. they release and receive theinformation particles perfectly), and that the information particles are not lost in the environment, the only factor that effects the PMF f Y ( T ) | X ( y ( T ) | x ) is the random motion of the MTs. Althoughthis assumption may not hold in practice, by assuming perfect transmitter and receiver, we onlyfocus on the effects of channel-shape on the information rate. Moreover, for some applicationssuch as pure mass transport in lab-on-chip devices, we are only concerned with mass transport asopposed to a communication system. The PMF f Y ( T ) | X ( y ( T ) | x ) could be used as an importantperformance measure in both cases. For a given TPCU duration T , as the number of MT tripsbetween the transmitter and the receiver increases, the number of particles that could potentiallybe transported increases. This results in higher mass transfer and potentially higher achievableinformation rates. Based on this observation we develop an optimization model which can beused to find the optimal channel shape.IV. C HANNEL S HAPE O PTIMIZATION M ODEL
Let the set G be the set of all possible cross sectional geometric shapes the channel couldhave, where the cross sectional shape is the shape of the channel when viewed from the top.In this paper we fist assume this set contains all the rectangular, regular polygonal, and regularpolygonal ring shapes and then extend the result to all two dimensional shapes. Our goal is tofind the optimal shape in these sets, given system parameters such as TPCU, concentration ofMTs, and average speed of the MTs. To setup this optimization problem, we model the effectsof channel shape on the number of MT trips.As discussed in the previous section, the number of MT trips during the TPCU interval has adirect effect on transportation of information particles. Without loss of generality let g ∈ G bethe channel under consideration. Therefore, let K ( T ) ( g ) be the number of trips during the TPCUinterval T inside the channel shape g , where a single MT trip is defined as the movement ofthe MT from anywhere in the channel to the transmission zone and then the receiver zone. Forexample, a single MT trip is shown in Figure 3. After the MT completes its first trip, subsequenttrips are defined as the movement of the MT from the receiver zone to the transmission zone andback. During any trip, an MT can deliver 0 or more information particles (up to its maximumload capacity). Because the motion of the MT filaments are random in nature, K ( T ) ( g ) is arandom variable. Therefore, to derive our optimization model we use the expected value of K ( T ) ( g ) .Let the random variable K ( T ) s ( g ) be the number of trips for a single MT during the TPCU
500 1000 1500 2000123456789 A v e r age N u m be r o f T r i p s TPCU (s) Poly, r=20, n=8, SimPoly, r=20, n=8, EstRing, ri=17, ro=22, n=20, SimRing, ri=17, ro=22, n=20, EstRec, w=20, l=60, SimRec, w=20, l=60, Est
Fig. 4: Average number of trips approximation compared with Monte Carlo simulations.interval, T . Let v avg be the average speed of the MTs, P ( g ) be the perimeter of the channelshape g ∈ G . A good estimate for the average number of MT trips, when a single MT is insidethe channel is given by E [ K ( T ) s ( g )] ≈ v avg TP ( g ) . (8)This approximation is based on the observation that on average from each trip, the MT travelsa distance equivalent to the perimeter of the channel. This assumption is verified using ourMonte Carlo simulation environment. The results of the simulations are shown in Figure 4. Inthis figure, we have considered a sample shape from each of our three different shape classes,and have shown that the approximations in Equation (8) (solid lines) are a good estimate of theactual average number of MT trips (point plots).In practice, there are typically more than one MT inside the channel. Furthermore, the numberof MTs is not a constant, but is dependent on the volume of fluid in the channel: the DNA-covered MTs are prepared in chemical solutions, and therefore have a constant concentrationinside the solution. Let C be the concentration of MTs as number of MTs per unit volume,and h be the height of the channel, and A ( g ) be the cross sectional area of the channel with geometric shape g ∈ G . Then the number of MTs inside the channel is given by M ( g ) = A ( g ) × h × C, (9)where M ( g ) is the number of MTs inside channel g .For each MT in the channel let the random variables K ( T ) i ( g ) be the number of trips during T seconds inside the channel with shape g , for i = 1 , , · · · , M . Then, the total number of MTtrips during T seconds is given by K ( T ) ( g ) = M (cid:88) i =1 K ( T ) i ( g ) . (10)The average number of MT trips during T seconds is therefore calculated as E [ K ( T ) ( g )] = M (cid:88) i =1 E [ K ( T ) i ( g )] . (11)We assume that the number of trips for individual MTs are independent and identically dis-tributed, because they are chemically similar and don’t interact with each other. Generally, thisassumption has been found to be valid in laboratory experiments [16], [30]. Due to the identicaldistributions, the equation simplifies to E [ K ( T ) ( g )] = M (cid:88) i =1 E [ K ( T ) s ( g )] = M ( g ) E [ K ( T ) s ( g )] . (12)Using the approximation shown in Equation (8), and Equation (9), the total number of trips canbe estimated as E [ K ( T ) ( g )] ≈ T × v avg × C × h × A ( g ) P ( g ) , (13)where A ( g ) and P ( g ) are the cross sectional area and perimeter of the channel with geometricshape g ∈ G , respectively. The channel shape optimization problem can then be formulated as max g ∈G E [ K ( T ) ( g )] . (14)Assuming all the other parameters are constant except the cross sectional shape (including theheight of the channel), the optimization problem becomes max g ∈G (cid:20) A ( g ) P ( g ) (cid:21) . (15)Equation (15) is very interesting because it states that the optimal shape is the one with thelargest area to perimeter ratio. However, there is an important constraint that must be satisfied. The perimeter of the channel must be small enough such that a single MT can complete at leasta single trip. Therefore, we assume the perimeter must be such that its length can be traveledby the MT, on average, during the TPCU interval T . This constraint ensures that the perimeteris small enough such that on average enough information particles can be delivered during thegiven TPCU duration. Finally, the optimization problem can be written in its complete form as max g ∈G (cid:20) A ( g ) P ( g ) (cid:21) subject to P ( g ) ≤ T v avg . (16)V. O PTIMAL S HAPE A NALYSIS
Using the optimization formula derived in the previous section, we analyze the optimal channeldesign for each shape class individually, and then the optimal design over all.
1) Rectangular Channels:
Let G Rec be the set of all rectangular shapes. Because any rectan-gular shape can be characterized by the two parameters width w and length l , the Equation (16)becomes max ( w,l ) ∈G Rec (cid:20) w × l w + 2 l (cid:21) subject to w + 2 l ≤ T v avg . (17)Solving this optimization problem, the optimal channel design given both TPCU interval T andaverage speed of the MTs v avg , is w = l = 0 . T v avg . Therefore, for rectangular channels, theoptimal channel shape is always the square shaped channel.
2) Regular Polygonal Channels:
Let G Poly be the set of all regular polygons. Regular polygonscan be characterized by two parameters: the number of sides n , and the radius of the circum-scribed circle r . Note that the set of all square shapes are also contained within the set G Poly (i.e. when n = 4 ). Using these parameters Equation (16) becomes max ( n,r ) ∈G Poly (cid:20) . nr sin(2 π/n )2 nr sin( π/n ) (cid:21) subject to nr sin( π/n ) ≤ T v avg . (18) Simplifying this equation we get max ( n,r ) ∈G Poly . r cos( π/n ) subject to nr sin( π/n ) ≤ T v avg , (19)where we have used the fact that sin(2 u ) = 2 sin( u ) cos( u ) . Based on this equation the optimalchannel is the one with n = ∞ and r = T v avg / π . Therefore, when regular polygons areconsidered the optimal shape is the circular-shaped channels. Because squares are also in the setof all regular polygons, we conclude that circular-shaped channels are also better than square-shaped channels.
3) Regular Polygonal Ring Channels:
Let G Ring be the set of all polygonal ring shapes. Theelements of this set can be characterized using three parameters: the number sides n , the radius ofthe inner polygon’s circumscribed circle r i , and the radius of the outer polygon’s circumscribedcircle r o . Note that G Ring ⊃ G
Poly , which means that the set G Ring contains all the regular polygonalshapes (this is achieved by setting r i = 0 ). Using these parameters the optimization problembecomes max ( n,r i ,r o ) ∈G Ring (cid:20) . n ( r o − r i ) sin(2 π/n )2 nr o sin( π/n ) (cid:21) subject to nr o sin( π/n ) ≤ T v avg . (20)Simplifying this equation we get max ( n,r i ,r o ) ∈G Ring (cid:20) . (cid:18) r o − r i r o (cid:19) cos( π/n ) (cid:21) subject to nr o sin( π/n ) ≤ T v avg . (21)Solving Equation 21, the optimal channel has parameters n = ∞ , r o = T v avg / π , and r i = 0 .This means the optimal channel is the circular-shaped channel. This may seem surprising atfirst, since one would expect that a ring-shaped circular channel would be better than a solidcircular-shaped channel at transporting information particles. However, because typically MTsfollow the walls of the channel, and that the number of MTs in a channel is proportional to thevolume of the channel [32], it becomes apparent that the solid circular channel would be betterthan a ring-shaped channel.
4) Overall Optimal Channel Design:
From the solution to the optimization formulas for eachof the three different shape classes, it is apparent that the optimal channel shape is the circular-shaped channel when the set of geometric shapes G RPR = G Rec ∪ G
Poly ∪ G
Ring is considered.Moreover, based on Equation (16) the optimal channel shape must have the largest area toperimeter ratio. Therefore, if the perimeter of the channel is fixed such that the constraint in theoptimization is satisfied with an equality, the optimal channel shape would be a circle for alltwo dimensional geometric shapes.VI. R
ESULTS AND D ISCUSSIONS
To verify our optimization formulas and their results, we rely on Monte Carlo simulation. Forthese simulations we use the same parameters discussed in Section II-B. Moreover, we assumethe height of the channel is always a constant h = 10 µ m regardless of the cross sectional shapeof the channel. The concentration of the MTs is also assumed to be C = 0 . MT/fL unlessspecified otherwise. The number of MT in the channel is always assumed to be
M T = (cid:98) A ( g ) hC (cid:99) , (22)unless specified otherwise.For the performance measure we rely on channel capacity [37], which is the maximum rateat which information can be transmitted reliably. Channel capacity is given by C = max f X ( x ) I ( X ; Y ( T ) ) , (23)where I ( X ; Y ( T ) ) is the mutual information between X and Y ( T ) . Mutual information is definedas I ( X ; Y ( T ) ) = E (cid:34) log f Y ( T ) | X ( y ( T ) | x ) (cid:80) x f Y ( T ) | X ( y ( T ) | x ) f X ( x ) (cid:35) , (24)where f Y ( T ) | X ( y ( T ) | x ) represents the probability of receiving symbol y ( T ) at the destinationduring TPCU interval T , given that symbol x was transmitted by the source; f X ( x ) representsthe probability of transmitting symbol x ; and E [ · ] represents expectation. Using Monte Carlosimulations the PMF f Y ( T ) | X ( y ( T ) | x ) is estimated for a given channel shape. For more detailsabout the Monte Carlo simulations for estimating this PMF please see [9], and for analyticalestimation of this PMF using Markov chain models please refer to [38]. In this work, we useMonte Carlo simulations since they are more accurate in multi-MT environments. Using the
20 25 30 35 40 45 50 20 25 30 35 40 45 500.80.911.11.21.31.41.51.61.71.8 Channel Width ( μ m)TPCU = 160 s, x max = 40Channel Length ( μ m) C hanne l C apa c i t y ( b i t s ) (a)
20 25 30 35 40 45 50 20 25 30 35 40 45 501.71.81.922.12.22.32.42.5 Channel Width ( μ m)TPCU = 240 s, x max = 40Channel Length ( μ m) C hanne l C apa c i t y ( b i t s ) (b) Fig. 5: Channel capacities of different rectangular channel shapes. The channel with the highestcapacity is shown using the blue dot. For TPCU of 160 seconds the optimal channel is 20 µ m × µ m (a), and for TPCU value of 240 seconds the optimal channel dimensions is 30 µ m × µ m (b).estimated PMF f Y ( T ) | X ( y ( T ) | x ) , the maximizing PMF f X ( x ) and hence the channel capacity iscalculated using Blahut-Arimoto algorithm [39], [40].We now consider each of three shape classes discussed in the previous sections, and use theEquations (17), (19), and (21) to find the optimal channel dimensions for each case. These resultsare then compared to the channel capacities obtained from Monte Carlo simulations. A. Rectangular Channels
First, we consider rectangular channels. We consider two values of TPCU: T = 160 seconds,and T = 240 seconds. Using the solution to Equation (17), the optimal channels are squarechannels with dimensions w = l = 0 . × × . µ m and w = l = 0 . × × . µ m for the T = 160 and T = 240 , respectively. To verify this result, we simulate all therectangular channels with widths and lengths raging from µ m to µ m in µ m steps.We assume x max (i.e. the maximum number of particles the transmitter can release) is 40. The C hanne l C apa c i t y ( b i t s ) x max (maximum number of possible transmission particles)TPCU value of 320 seconds Poly n = 20, r = 25.57, MT = 20Poly n = 16, r = 25.63, MT = 16Poly n = 12, r = 25.76, MT = 19Poly n = 8, r = 26.13, MT = 19Rec w = 40, l = 40, MT = 16 Fig. 6: Channel capacity of regular polygonal channels with constant perimeter of 160 µ mversus maximum number of particles that can be released at the transmitter.channel capacity is then calculated based on the simulations. Figure 5 shows the results for T = 160 seconds (a), and T = 240 seconds (b). As can be seen the channels with the largestchannel capacity are the 20 µ m × µ m, and 30 µ m × µ m, respectively. These results arein perfect agreement with the results obtained from our optimization model. B. Regular Polygonal Channels
In Section V we showed that for regular polygonal channels the optimal shape is a circular-shaped channel. To verify this result, we first consider a constant TPCU value of 320 seconds.From Equation (16), the maximum channel perimeter that satisfies the constraint is µ m.Therefore, we assume that the channel perimeter is fixed at this value and show that the capacityincreases a the channel shape become more circular by increasing the number of sides n . Inparticular we consider the following channels: square channel of length 40 µ m (equivalent toa polygon with n = 4 ), and regular polygonal channels with parameters ( n = 8 , r = 26 . µ m), ( n = 12 , r = 25 . µ m), ( n = 16 , r = 25 . µ m), and ( n = 20 , r = 25 . µ m). C hanne l C apa c i t y ( b i t s ) x max (maximum number of possible transmission particles)TPCU value of 250 seconds Poly n = 20, r = 17, MT = 8Poly n = 20, r = 20, MT = 12Poly n = 20, r = 22.75, MT = 16Poly n = 20, r = 25.57, MT = 20 Fig. 7: Channel capacity of four circular-shaped channels versus maximum number of particlesthat can be released at the transmitter.Based on our optimization formula we expect the 20-sided channel to be the optimal among allthese channels since it is more circular. Figure 6 shows the channel capacity for each of thesechannels obtained through Monte Carlo simulations versus the maximum number of particles thetransmitter can release x max . As can be seen, the optimal channel is indeed the 20-sided regularpolygon. Moreover, from the obtained pattern it is obvious that as the number of sides increasethe channel capacity increases. This supports the results obtained from our optimization formulathat circular channels are optimal.To evaluate the performance of Equation (19) we consider the TPCU interval of T = 250 seconds. The solution of this optimization problem is a channel with parameters n = ∞ and r ≈ µ m. Because in our simulation environment the number of sides must be a finite number,we use the value of n = 20 to simulate a perfect circle. We then simulate four different regularpolygonal channels with parameters: ( n = 20 , r = 17 µ m), ( n = 20 , r = 20 µ m), ( n = 20 , r = 22 . µ m), and ( n = 20 , r = 25 . µ m). Figure 7 shows the channel capacity of eachchannel versus the maximum number of particles the transmitter could release. As can be seen C hanne l C apa c i t y ( b i t s ) x max (maximum number of possible transmission particles)TPCU value of 320 seconds Ring n = 20, ro = 25.57, ri = 0, MT = 20Ring n = 20, ro = 25.57, ri = 10, MT = 17Ring n = 20, ro = 25.57, ri = 15, MT = 13Ring n = 20, ro = 25.57, ri = 20, MT = 7 Fig. 8: Channel capacity of ring-shaped channels versus maximum number of particles that canbe released at the transmitter.in the figure, the channel with the radius of 20 µ m achieves the highest capacity among the fourchannels. This result is in agreement with the optimal solution derived from our optimizationformula. C. Regular Polygonal Ring Channels
In this section, we validate our optimization formula for polygonal ring-shaped channelspresented in Equation (21) through simulations. We assume the TPCU interval is fixed at 320seconds. Using the solution to Equation (21), the optimal channel dimensions are n = ∞ , r i = 0 , and r o ≈ . µ m. In Figure 6 we have already shown that the solid circular-shapedchannel is optimal among regular polygons. Therefore, we only consider ring-shaped channelsfor validation. In particular, we consider the following ring-shaped channels: ( n = 20 , r o = 25 . µ m, r i = 0 µ m), ( n = 20 , r o = 25 . µ m, r i = 10 µ m), ( n = 20 , r o = 25 . µ m, r i = 15 µ m), and ( n = 20 , r o = 25 . µ m, r i = 20 µ m). Figure 8 shows the channel capacity of thesechannels versus the maximum number of particles the transmitter could release. We can see from the figure that a solid circular channel achieves the highest capacity as predicted by ouroptimization formula.Based on all these results, we can see that a solid circular-shaped channel is the optimal shapethat maximizes information rate, when stationary kinesin with mobile MT filaments are used fortransportation in molecular communication systems. Moreover, as the number of sides increaseand the channel becomes closer to circular, the performance gains of adding more sides decreases.Based on the obtained results, there is negligible performance difference between channels withthe number of sides greater than 20. Therefore, these channels are sufficient representation ofperfectly circular-shaped channels. VII. C ONCLUSION
In this work we considered the problem of finding the optimal channel shape for on-chipmolecular communication based on kinesin driven microtubule (MT) motility. To setup ouroptimization model, we first derived a mathematical model relating the average number of MTtrips to the shape of the channel. Using this model we then presented the three optimizationformulas and their respective solutions for three different classes of shapes: the rectangular,regular polygonal, and regular polygonal rings. We showed that the optimal solution is thesquare channel for the rectangular channel shapes, and circular shape for the regular polygonaland ring-shaped channels. Finally, we showed that when all classes of 2 dimensional shaped areconsidered the circular-shaped channel tend to be the one that achieves the highest informationrate. Moreover, using our optimization formula we showed that the dimensions of the optimalchannel can be reliably estimated. R
EFERENCES [1] R. Daw and J. Finkelstein, “Insight: Lab on a chip,”
Nature , vol. 442, no. 7101, pp. 367–418, 2006.[2] P. Yager, T. Edwards, E. Fu, K. Helton, K. Nelson, M. R. Tam, and B. H. Weigl, “Microfluidic diagnostic technologiesfor global public health,”
Nature , vol. 442, pp. 412–418, July 2006.[3] I. F. Akyildiz, F. Brunetti, and C. Blazquez, “Nanonetworks: A new communication paradigm,”
Computer Networks ,vol. 52, no. 12, pp. 2260–2279, 2008.[4] S. F. Bush,
Nanoscale Communication Networks . Artech House, 1 ed., 2010.[5] S. Hiyama, Y. Moritani, T. Suda, R. Egashira, A. Enomoto, M. Moore, and T. Nakano, “Molecular communication,” in
Proc. of 2005 NSTI Nanotechnology Conference , pp. 391–394, 2005.[6] S. Hiyama and Y. Moritani, “Molecular communication: harnessing biochemical materials to engineer biomimeticcommunication systems,”
Nano Communication Networks , vol. 1, no. 1, pp. 20–30, 2010.[7] A. W. Eckford, “Timing information rates for active transport molecular communication,” in
Nano-Net , vol. 20 of
LectureNotes of the Institute for Computer Sciences, Social Informatics and Telecommunications Engineering , pp. 24–28, SpringerBerlin Heidelberg, 2009.[8] N. Farsad, A. W. Eckford, S. Hiyama, and Y. Moritani, “Quick system design of vesicle-based active transport molecularcommunication by using a simple transport model,”
Nano Communication Networks , vol. 2, no. 4, pp. 175–188, 2011. [9] N. Farsad, A. Eckford, S. Hiyama, and Y. Moritani, “On-Chip Molecular Communication: Analysis and Design,” IEEETransactions on NanoBioscience , vol. 11, no. 3, pp. 304–314, 2012.[10] M. U. Mahfuz, D. Makrakis, and H. T. Mouftah, “On the characterization of binary concentration-encoded molecularcommunication in nanonetworks,”
Nano Communication Networks , vol. 1, no. 4, pp. 289–300, 2010.[11] L. C. Cobo and I. F. Akyildiz, “Bacteria-based communication in nanonetworks,”
Nano Communication Networks , vol. 1,no. 4, pp. 244–256, 2010.[12] T. Nakano, T. Suda, T. Koujin, T. Haraguchi, and Y. Hiraoka, “Molecular communication through gap junction channels,”
Springer Transactions on Computational Systems Biology X , vol. 5410, pp. 81–99, 2008.[13] M. Pierobon and I. F. Akyildiz, “A physical end-to-end model for molecular communication in nanonetworks,”
IEEEJournal on Selected Areas in Communications , vol. 28, no. 4, pp. 602–611, 2010.[14] M. Pierobon and I. Akyildiz, “Capacity of a diffusion-based molecular communication system with channel memory andmolecular noise,”
IEEE Transactions on Information Theory , vol. 59, no. 2, pp. 942–954, 2013.[15] K. V. Srinivas, A. Eckford, and R. Adve, “Molecular communication in fluid media: The additive inverse gaussian noisechannel,”
IEEE Transactions on Information Theory , vol. 58, no. 7, pp. 4678–4692, 2012.[16] S. Hiyama, Y. Moritani, R. Gojo, S. Takeuchi, and K. Sutoh, “Biomolecular-motor-based autonomous delivery of lipidvesicles as nano- or microscale reactors on a chip,”
Lab on a Chip , vol. 10, no. 20, pp. 2741–2748, 2010.[17] M. Gregori and I. Akyildiz, “A new nanonetwork architecture using flagellated bacteria and catalytic nanomotors,”
IEEEJournal on Selected Areas in Communications , vol. 28, no. 4, pp. 612–619, 2010.[18] P. Lio and S. Balasubramaniam, “Opportunistic routing through conjugation in bacteria communication nanonetwork,”
Nano Communication Networks , vol. 3, pp. 36–45, Mar. 2012.[19] M. J. Moore, T. Suda, and K. Oiwa, “Molecular communication: modeling noise effects on information rate,”
IEEETransactions on NanoBioscience , vol. 8, no. 2, pp. 169–180, 2009.[20] A. Enomoto, M. J. Moore, T. Suda, and K. Oiwa, “Design of self-organizing microtubule networks for molecularcommunication,”
Nano Communication Networks , vol. 2, no. 1, pp. 16–24, 2011.[21] T. Nakano, M. Moore, F. Wei, A. Vasilakos, and J. Shuai, “Molecular communication and networking: Opportunities andchallenges,”
IEEE Transactions on NanoBioscience , vol. 11, pp. 135–148, june 2012.[22] M. Briani, G. Germani, E. Iannone, M. Moroni, and R. Natalini, “Design and optimization of reaction chamber anddetection system in dynamic labs-on-chip for proteins detection,”
IEEE Transactions on Biomedical Engineering , vol. 60,pp. 2161–2166, Aug 2013.[23] S. Alla, A. Huddle, J. Butler, P. S. Bowman, J. Clark, and F. Beyette, “Point-of-care device for quantification of bilirubinin skin tissue,”
IEEE Transactions on Biomedical Engineering , vol. 58, pp. 777–780, March 2011.[24] Y. Chahibi, M. Pierobon, S. O. Song, and I. Akyildiz, “A molecular communication system model for particulate drugdelivery systems,”
IEEE Transactions on Biomedical Engineering , vol. 60, pp. 3468–3483, Dec 2013.[25] T. Kim, M.-T. Kao, E. F. Hasselbrink, and E. Meyhfer, “Active alignment of microtubules with electric fields,”
NanoLetters , vol. 7, no. 1, pp. 211–217, 2007.[26] I. Dujovne, M. van den Heuvel, Y. Shen, M. de Graaff, and C. Dekker, “Velocity modulation of microtubules in electricfields,”
Nano Letters , vol. 8, no. 12, pp. 4217–4220, 2008.[27] E. Kim, K.-E. Byun, D. S. Choi, D. J. Lee, D. H. Cho, B. Y. Lee, H. Yang, J. Heo, H.-J. Chung, S. Seo, and S. Hong,“Electrical control of kinesinmicrotubule motility using a transparent functionalized-graphene substrate,”
Nanotechnology ,vol. 24, p. 195102, May 2013.[28] S.-Y. Teh, R. Lin, L.-H. Hung, and A. P. Lee, “Droplet microfluidics,”
Lab Chip , vol. 8, pp. 198–220, 2008.[29] T. Nitta, A. Tanahashi, M. Hirano, and H. Hess, “Simulating molecular shuttle movements: Towards computer-aided designof nanoscale transport systems,”
Lab on a Chip , vol. 6, no. 7, pp. 881–885, 2006.[30] T. Nitta, A. Tanahashi, and M. Hirano, “In silico design and testing of guiding tracks for molecular shuttles powered bykinesin motors,”
Lab on a Chip , vol. 10, no. 11, pp. 1447–1453, 2010.[31] N. Farsad, A. W. Eckford, S. Hiyama, and Y. Moritani, “A simple mathematical model for information rate of activetransport molecular communication,” in
Proc. of 2011 IEEE INFOCOM Workshops , (Shanghai, P. R. China), pp. 473–478,2011.[32] J. Clemmens, H. Hess, J. Howard, and V. Vogel, “Analysis of microtubule guidance in open microfabricated channelscoated with the motor protein kinesin,”
Langmuir , vol. 19, no. 5, pp. 1738–1744, 2003.[33] N. Farsad, A. W. Eckford, and S. Hiyama, “A mathematical channel optimization formula for active transport molecularcommunication,” in
Proc. of 2012 IEEE International Conference on Communications (ICC) Workshops , (Ottawa, Canada),2012.[34] N. Farsad, A. Eckford, and S. Hiyama, “Modelling and design of polygon-shaped kinesin substrates for molecularcommunication,” in
Proc. of 12th IEEE Conference on Nanotechnology (IEEE-NANO) , (Birmingham, UK), pp. 1–5, 2012.[35] D. Steuerwald, S. M. Fruh, R. Griss, R. D. Lovchik, and V. Vogel, “Nanoshuttles propelled by motor proteins sequentiallyassemble molecular cargo in a microfluidic device,”
Lab Chip , vol. 14, pp. 3729–3738, 2014.[36] S. Hiyama, R. Gojo, T. Shima, S. Takeuchi, and K. Sutoh, “Biomolecular-motor-based nano- or microscale particletranslocations on DNA microarrays,”
Nano Letters , vol. 9, no. 6, pp. 2407–2413, 2009.[37] T. M. Cover and J. A. Thomas,
Elements of Information Theory 2nd Edition . Wiley-Interscience, 2 ed., 2006. [38] N. Farsad, A. Eckford, and S. Hiyama, “A markov chain channel model for active transport molecular communication,” IEEE Transactions on Signal Processing, , vol. 62, pp. 2424–2436, May 2014.[39] R. Blahut, “Computation of channel capacity and rate-distortion functions,”
IEEE Transactions on Information Theory ,vol. 18, no. 4, pp. 460–473, 1972.[40] S. Arimoto, “An algorithm for computing the capacity of arbitrary discrete memoryless channels,”