A molecular communications model for drug delivery
aa r X i v : . [ c s . ET ] A ug A molecular communications model for drugdelivery
Mauro Femminella, Gianluca Reali, Athanasios V. Vasilakos ∗ Abstract — This paper considers the scenario of a targeteddrug delivery system, which consists of deploying a numberof biological nanomachines close to a biological target (e.g. atumor), able to deliver drug molecules in the diseased area.Suitably located transmitters are designed to release a continuousflow of drug molecules in the surrounding environment, wherethey di ff use and reach the target. These molecules are receivedwhen they chemically react with compliant receptors deployedon the receiver surface. In these conditions, if the release rateis relatively high and the drug absorption time is significant,congestion may happen, essentially at the receiver site. Thisphenomenon limits the drug absorption rate and makes thesignal transmission ine ff ective, with an undesired di ff usion ofdrug molecules elsewhere in the body. The original contributionof this paper consists of a theoretical analysis of the causes ofcongestion in di ff usion-based molecular communications. For thispurpose, it is proposed a reception model consisting of a setof pure loss queuing systems. The proposed model exhibits anexcellent agreement with the results of a simulation campaignmade by using the Biological and Nano-Scale communicationsimulator version 2 (BiNS2), a well-known simulator for molecu-lar communications, whose reliability has been assessed throughin-vitro experiments. The obtained results can be used in ratecontrol algorithms to optimally determine the optimal releaserate of molecules in drug delivery applications. Index Terms — Drug delivery, molecular communications, dif-fusion, congestion, service time, queuing model.
I. I ntroduction
Molecular communication is a novel paradigm allowinginformation exchange between biological nanomachines (bio-nanomachines) over short ranges, within an aqueous envi-ronment [1]. In this model, a transmitter nanomachine emitsmolecules, which propagate and eventually are received bya receiver nanomachine. The information is usually encodedin either the timing or the amplitude of the concentration ofmolecules. The reception process usually consists of a chem-ical reaction between those molecules (ligand) and compliantreceptors present on the receiver surface. Bio-nanomachinesare made of biological materials and can perform a numberof simple tasks which need to be coordinated, by means ofmessage exchange, in order to accomplish complex functions.Di ff usion-based molecular communications are a specialform of molecular communications [2], very common in This work was supported in part by EU project H2020 FET Open CIRCLE(Coordinating European Research on Molecular Communications). Asteriskindicates corresponding author.M. Femminella, and G. Reali are with the Department of Engineering,University of Perugia, via G. Duranti 93, 06125 Perugia, Italy, (email:{mauro.femminella,gianluca.reali}@unipg.it).A. V. Vasilakos is with Department of Computer Science, Electricaland Space Engineering, Lulea University of Technology, Sweden (email:[email protected]). nature. They can be roughly classified in two categories,namely pure-di ff usion and di ff usion with drift. For instance,the former includes communications in the extracellular matrixof the connective tissue, whilst the latter encompasses com-munications in circulatory and lymphatic systems.In this paper, we consider drug delivery in pure-di ff usionmolecular communications. In drug delivery systems a nearlycontinuous flow of molecules is transmitted from the transmit-ter (TX) to the receiver (RX) [3], in order to achieve a desirede ff ective average drug delivery rate at a specific target. In theseapplications, the TX node has the role of an actuator, and theRX node is either the target or a sensor which measures andcontrols the delivery of drug / grow factors molecules to thesurrounding cells.In this scenario, the space around the receiver could becomefull of drug molecules, which produce congestion. Congestioncontrol in packet networks is a research topic that attractedresearchers in the last 30 years (see, e.g., [4], [5]), since theseminal works of Van Jacobson [6] and K.K. Ramakrishnan[7]. However, until now, little attention has been paid tothe issue of congestion in molecular communications, be-yond some initial investigations reported in [8], [9]. Indeed,molecular communications, especially those making use ofdi ff usion-based propagation [10], are not only very sensitiveto congestion, but also di ffi cult to control due to the highcommunication latency. Thus, an early detection of congestionis necessary in order to avoid heavy receiver overload. How-ever, beforehand, we have to give a definition of congestionin molecular communications.In molecular communications, congestion consists in theinability of receivers to capture all molecules in the receptorspace. This limitation is due to the combination of two factors,which are the number of receptors present over the receiversurface and the tra ffi cking time [11], which can be regardedas the molecule reception time. The tra ffi cking time is definedas the time needed to make a binding between a ligand anda compliant receptor, plus the time to internalize the resultingcomplex, plus the time to re-expose another receptor. It canlast up to tens of seconds [11] and significantly a ff ect themolecule absorption process. Motivation : in the bio-medical field, the receiver congestionproblem is referred to as “receptor saturation”. Building amodel for receptor saturation is a very important open issue,since the lack of a general model requires executing expensivelab experiments for receptors characterization [12]. In addition,receptors saturation a ff ects drug delivery and disease modeling[13]. However, most of research in the molecular commu-nications either neglect the presence of individual receptors(absorbing receiver [14]), or neglect the absorption time of molecules (absorbing receptors [15]). A more accurate model,which partially addresses these shortcomings by modelingmolecules absorption via receptors with a birth-death Markovprocess, is presented in [16]. However, also this model has itsshortcomings. In fact, it models receptors as a pool of servers,which can be invoked and used when any new molecule entersthe system. Instead, receptors are isolated service elements,and this nature must be taken into account. Main contributions : the paper has two main contributions.The first one is the analysis of the causes of congestion atRX when a continuous flow of molecules is transmitted bythe TX. The second and main contribution is a queuing modelthat includes these causes, which has been validated by anextensive simulation campaign. Our idea is to model eachreceptor on the RX surface as a pure loss queuing system.Specifically, we model each receptor by an M / M / / ffi cking time.Given the unguided di ff usion of the transmitted signal, wemodel the receiving nanomachine as a node with a largenumber of receiving antennas (receptors), each interfaced toa layer 2 with very high processing time (the tra ffi ckingtime) without waiting room. Thus, when the concentration ofmolecules nearby the RX surface becomes excessive (channelcongestion) and a significant part of receptors are busy, alsothe RX becomes congested (receiver congestion). This anal-ogy between di ff usion-based molecular communications andwireless communications is motivated also by recent papersdealing with multiple-input and multiple-output schemes [17]and directional antennas [18] to build directional receivers inthe field of molecular communications, similarly as in thewireless world.Finally, although the queuing theory has already been usedin the field of bio-physics for modeling ligand-receptor sys-tems [11], [19], the novel contribution of this paper consistsof using this theory for modeling receptor saturation from thecommunications perspective, thus treating it as a congestionproblem in communications networks. To the best of ourknowledge, this is the first paper using this approach.The paper is organized as follows. In section II, we illustratethe related work in the field. Section III presents the systemmodel, including a detailed description of our proposed model,which is able to capture the root causes of congestion. Nu-merical results, which include the results of the simulationcampaign, used to validate the proposed mechanisms, arepresented in Section IV. Finally, we draw our conclusions inSection V. II. B ackground and R elated W orks A. Molecular communications
A complete view of a layered network architecture formolecular communications has been proposed in [20]. Follow-ing the layered architecture of traditional communication net-works, such as the Open Systems Interconnection model (OSI)and TCP / IP reference model, a formal model was developedfor each layer, and potential research directions are illustratedfor each of them. Other works explore issues related to thehigh protocol layers, including connection-oriented protocols [9], feedback-based rate control [8], network coding [21], andsecurity issues [22].However, major research e ff orts in molecular communica-tions are focused on the physical layer issues of di ff erentcommunication media. In particular, the information capacityand the physical features (e.g., delay, signal attenuation andamplification) of molecular communications are studied by us-ing random walk models [23]–[25], random walk models withdrift [26], [27], di ff usion-based models [28]–[32], di ff usion-reaction-based models [33], [34], active transport models [23],[35], and a collision-based model [36].As for the di ff usion-based models, which is the modelconsidered in this work, a review of di ff erent transmissionschemes is provided in [37], where the existing schemesare classified into pulse position modulation (PPM, [26]),concentration shift keying (CSK, [38]), and molecule shift key(MoSK), a further proposed category which combines burstsof di ff erent molecules to encode data. The possible trade-o ff between symbol duration and communication distance isanalyzed in [39].Some recent papers show enhanced receivers for di ff usion-based molecular communications [25], [39]–[45]. However,these works, as most of research papers in molecular commu-nications, do not take into account that signal reception occursvia binding between a signal molecule and the relevant surfacereceptor, and the number of receptors is typically assumed tobe so large to cover all the receiving nanomachine surface.This is referred as the so-called absorbing receiver model [14],in which a molecule, when gets in touch with the surfaceof the RX node, is immediately absorbed and removed fromthe surrounding environment. Although the absorbing receivermodel is acceptable for some processes used by cells to engulfother cells, microvesicles, or proteins (phagocytosis and / orpynocytosis, [46]), most of cells communications are carriedout via receptors. Thus, in some situations the assumption ofusing an absorbing receiver is a very rough approximation,since receptors are a finite resource and each nanomachineusually has a limited, yet large, number of receptors ofa given type. A variation of the absorbing receiver modelis the receiver with absorbing receptors [15], in which amolecule is instantaneously absorbed when it hits a compliantreceptors. Again, although this model is acceptable in somespecific situations (e.g. for modeling ions channels [47]), inmost of cases it does not, since the stochastic nature of theligand-receptor binding is not considered. Di ff erently, it iscontemplated in [16], and further investigated in [48], whichalso includes the tra ffi cking time [11] in the reception process. B. Drug delivery systems
One of the most popular topics in nanomedicine is thetargeted drug delivery, since it could be the basis for themodern medical therapeutics. Its main goal is to provide alocalized drug delivery only where medication is needed, thusavoiding to a ff ect other healthy parts of the body. This is dueto the very small size of the particles delivering medicine,which is in the order of nanometers (nanoparticles), whichallow them to di ff use into the bloodstream and across the vascular and interstitial barriers. Nanoparticles present ligandson their surface, in order to improve cell targeting. In fact,these ligands increase the chances of binding to the surfacereceptors of target cells [49].Many types of drug delivery, based on molecular com-munications, have been proposed in the technical literature,with both passive and active transport of molecules. For whatconcerns passive transport, the most popular approach is basedon particle di ff usion, with or without drift. In case the drugdelivery happens via the circulatory or lymphatic system, themost appropriate model is the di ff usion with drift, by alsoconsidering the e ff ect of collisions with blood cells [50], [51].In case the drug delivery happens outside blood vessels, e.g.in the extracellular matrix, the mere di ff usion-based modelis su ffi ciently accurate. In case of active transport, the mainalternatives are those based on bacteria.In order to design an e ffi cient and targeted drug delivery, itis essential to correctly determine the distribution of particlesover time. For this reason, the relevant research has producedseveral approaches, from statistical models [52] to analyticalones [53]. In particular, drugs delivered through bloodstreamcan reach any part of the body, by passing through thecomplex network of blood vessels (Particulate Drug DeliverySystem, PDDS). In [53], further expanded in [54], [55], theauthors propose a specific model of the cardiovascular system,by analyzing the peculiarities of di ff erent vessels, from thesmallest to the largest ones, thus obtaining a drug propagationnetwork model which takes into account even bifurcations andjunctions of vessels. The resulting model includes also the in-dividual features of the cardiovascular system by physiologicalparameters, such as the heartbeat rate profile. This model couldbe useful to guarantee that the drug concentration does notreach toxic levels in the body and it is essential for optimizingthe drug injection, by identifying the suitable injection pointsthat guarantee the best delivery profile for maximizing thetreatment benefits while minimizing the amount of drugs inthe other parts of the body. To this aim, it is necessary torealize a complete control of the drug release, from the initialrelease, to the suitable release rate.A di ff erent and even more innovative model considers thedelivery from a number of nanomachines (nano-actuator),implanted close to the target (e.g. a tumor), which bypasses theinjection through the cardiovascular system, thus minimizingthe side e ff ects on the healthy parts of the body. A theoreticalmodel on the local rate control between a nano-controllerand a nano-actuator has been presented in [8]. However, thisis a high abstraction model, which does not consider finegrained dynamics and does not explore the receiver congestionphenomenon. Subsequently, a related communication protocol,inspired by the TCP congestion control, has been proposedin [9]. This protocol can control the drug release process byusing feedback messages. A challenging aspect of the protocoloperation is the usage of feedback messages, which makes theprotocol mostly suitable for environments without drift.Drug delivery systems lay their foundations on the interac-tions of drug molecules with the receptors at target cells [56]–[58]. Two opposite main theories exist, namely the occupancytheory [59] and the rate theory [60]. The former assumes that the drug e ff ects are maximized when all receptors at the targetare occupied by the drug, whereas the rate theory is developedon the assumption that the excitation of a stimulant drugis proportional to the rate of the drug-receptor combination,rather than to the proportion of receptors occupied by thedrug. However, behavior of drugs can be quite di ff erent. Forinstance, it is not always necessary to reach an occupancy of100% of receptors to produce a full response on the target cell[61]. In addition, it was found that a good model for describingthe e ff ects of the drug-receptor interaction depends on the typeof drug [62]. In any case, our proposed receiver model can beused for analyzing both of them.III. T he S ystem M odel The considered drug delivery scenario consists in twofixed bio-nanomachines, staying at a distance d between theircenters, as depicted in Fig. 1. One of them is the transmittingnode and the other is the receiver. The communication happensby means of signal molecules released by TX, which propagateby di ff usion, modeled as Brownian motion [10]. The Brownianmotion is characterized by the di ff usion coe ffi cient, given by D = K b T πη r c , tx . K b is the Boltzmann constant, T is the temperatureexpressed in Kelvins, η is the viscosity of the medium, and r c , tx is the radius of the considered molecules. The shape of bothRX and TX is spherical with radius r RX and r T X , respectively.The number of receptors deployed on the RX surface is R RX ,and they are modeled as a circular area of radius r r , rx .The transmission is organized in bursts of molecules of size Q , spaced in time by a short period equal to ∆ t , so as toemulate a continuous emission of drug molecules with rate Q / ∆ t . The reception process at receptors is modeled as anexponential random variable with mean T tra f f . We assumethat the communication session can be set up and torn downwith a specific protocol, e.g. similar to the one shown in [9].However, since our results are una ff ected by the protocol used,the latter is beyond the scope of this paper.Let us assume that the center of TX is located at the originof the system of coordinates, and initially we neglect thedisturbing presence of the RX. We recall that the RX is notabsorbing, although able to capture molecules. We now solvethe di ff usion equation upon the transmission an impulse ofsize Q molecules at t =
0. The concentration of molecules c ( t , r ), in a given point of the space at a distance r from thecenter of system of coordinates and at a specific time t , canbe obtained by solving the di ff usion equation (Fick’s secondlaw of di ff usion [10]), which is ∂ c ( t , r ) ∂ t = D ∇ c ( t , r ) , (1)obtaining the system impulse response scaled by Q : h Q ( t , r ) = Q (4 π Dt ) e (cid:18) − r Dt (cid:19) . (2)Since the di ff usion equation is a linear equation [63], thesolution to (1) for a train of bursts of size Q spaced by a TX RX d r RX a k r ( a k ) r TX S k a k-1 b Fig. 1: System model. We assume there is a circular symmetryaround the axis d .period ∆ t ( P ⌊ t / ∆ t ⌋ i = Q δ ( t − i ∆ t )) and started at t = c ( t , r ) = ⌊ t / ∆ t ⌋ X i = h Q ( t − i ∆ t , r ) ≈ ∆ t Z t h Q ( τ, r ) d τ −−−−→ t →∞ Q ∆ t π Dr . (3)This is a significant result, although obtained with the simpli-fying assumption of neglecting the presence of the RX node.The arrival rate of molecules in the receiver space dependson the concentration close to the receiver surface [16], which,in turn, from (3), results to be inversely proportional to thedistance from the source of molecules, that is c ( r ) ∝ r .As mentioned above, most of the receiver models do nottake into account the presence of receptors on the receiversurface, and assume that the receiver can capture all compli-ant molecules in the surrounding space (absorbing receiver).Clearly, these models are unsuitable for detecting congestionsince it is due to the limited rate of creation of ligand-receptorsbonds, and depends on both the number of receptors R RX andon the tra ffi cking time T tra f f .In this regard, a more accurate model is the “reversiblefirst-order reaction” model presented in [16]. Although thismodel has been proposed for a di ff erent goal, it is instructiveto review it in the perspective of congestion in molecularcommunications. In addition, it allows further detailing thesystem model. A. Review of reversible first-order reaction model
According to [16], the status of the receiver at any time t can be represented by the value of a random variable ˜ n ( t ),which denotes the number of ligand-receptor bonds at time t . Thus, the receiver behavior is modeled by the randomprocess { [˜ n ( t )] , t ≥ } . In general, ˜ n ( t ) increases when newligand-receptor bonds are created, and decreases due to theinternalization of these bonds [11], that is ligand moleculesassimilation. In our model, we neglect the probability of thesecondary e ff ect of the rupture of the ligand-receptor bond be-fore its internalization, in order to simplify the comprehensionof the main underlying phenomena. However, if needed, it canbe easily introduced in the receiver model.The random process { [˜ n ( t )] , t ≥ } is modeled as a birth-death, discrete-valued, continuous-time Markov process, with˜ n ( t ) ∈ { , , , . . ., R RX } [64]. Let π i ( t ) denote the probabilitythat the Markov chain is in state i at time t ; that is, π i ( t ) = P { ˜ n ( t ) = i } . This means that the system has i active ligand-receptor bonds. Then the dynamical equations for the state probabilities are as follows: ddt π i ( t ) = λ i − π i − ( t ) − ( λ i + µ i ) π i ( t ) + µ i + π i + ( t ) , (4)where is noted that all the rates λ i and µ i are state-dependent.In particular, we observe that • the coe ffi cient λ i describes the bond formation rate whenthe number of active bonds are i . In a queuing system, itrepresents the users arrival rate when the system is in state i . It depends upon the number of available receptors onthe cells surface, R RX − i , the concentration of moleculesclose to the cells surface at a given time, c ( t , r ), and onthe per receptor association rate constant, k + [16], [65].Specifically, the constant k + in [16] is defined as theproduct of Z , which represents how frequently a collisionbetween a particle and an unbound receptor happens,and F C ( E a ), which is the fraction of collisions havingan energy higher than the threshold E a , so that they canform a bond. Clearly, it results that Z is proportional to theaverage number of free receptors, that is Z ∝ ( R RX − n b ),where n b is the average number of bond receptors insteady-state conditions. In turn, n b is a function of theaverage arrival rate, λ a , formally defined in (16). As for F C ( E a ), in our model, each collision between a ligandand a compliant receptor corresponds to an assimilation,thus F C ( E a ) =
1. The resulting rate is the product of thesethree quantities, that is, λ i = k + c ( t , r ) ( R RX − i ) = Zc ( t , r ) ( R RX − i ) ; (5) • the coe ffi cient µ i describes the bond internalization ratewhen the number of active bonds are i . A memorylessmodel for this phenomenon is widely accepted in litera-ture [11]. It depends upon the number of currently formedcomplexes and on the average tra ffi cking time, wherethe tra ffi cking time is a random variable exponentiallydistributed with average value T tra f f . This means that µ i = µ i = i / T tra f f . (6)We assume that the TX node sends molecules according tothe constant pattern with rate Q / ∆ t . Then, in the steady state,the arrival process at the receiver can be approximated asa Poisson process, due to the randomness and memorylessproperties of Brownian motion. The average arrival rate λ a observed by the receiver is given by λ a = R RX − X i = λ i π i . (7)However, in addition to the arrival rate, also the value of therejection rate λ r is very important. In fact, since the modelconsists of a pure loss queuing system with R RX servers andstate dependent arrival rate, it clearly exhibits some losses. λ r represents the rate at which a ligand-receptor bond cannotbe formed, due to hits of ligand molecules with already busyreceptors. By following a reasoning similar to the one reportedabove for evaluating λ i , we can say that the state dependentrejection rate is equal to λ i , r = Z ∗ c ( t , r ) i = λ i n b i ( R RX − n b ) ( R RX − i ) . (8) Clearly, the per state rejection rate is proportional to thenumber of busy receptor, i , and to the Z ∗ , which represents howfrequently a collision between a particle and a busy receptorhappens, which, in turn, is obviously proportional to n b , andthus Z ∗ = Z n b R RX − n b . Consequently, the average rejection rate λ r can be evaluated as λ r = R RX X i = λ i , r π i . (9) B. Proposed congestion model
We model the reception process of any single receptor asan M / M / / α as r ( α ) = q ( d + R RX sin ( α )) + ( R RX cos ( α )) , (10)with α ∈ [ − π/ , π/ λ ∗ a ), which is a measurable parameter, it is thesuperposition of the rates of arrival to each receptor j , denotedas λ ∗ a , j . Each receptor j is identified by the corresponding angle α j and its distance from the TX center is equal to r ( α j ). Wedenote as j = α = π/
2, for whichthe distance to the center of TX is d + r RX . Considering that,as in the previously illustrated model [16], the arrival rateis proportional to the local concentration of molecules, andtaking into account (3), it results that λ ∗ a = R RX X j = λ ∗ a , j = λ ∗ a , ( d + r RX ) R RX X j = r ( α j ) , (11)where λ ∗ a is a system parameter which can be easily measuredby the receiver. Thus, the assimilation rate for each receptor λ ∗ a , j can be derived by λ ∗ a by inverting (11) and finding λ ∗ a , .In an M / M / / ffi c λ a ,the rate of the rejected tra ffi c is λ r = T traf f λ a − T traf f λ a . Thus, by usingthis result for each receptor, it is easy to find the system levelrejection rate λ ∗ r , given by λ ∗ r = R RX X j = λ ∗ r , j = λ ∗ a , ( d + r RX ) R RX X j = r ( α j ) (cid:18) r ( α j ) λ ∗ a , ( d + r RX ) T traf f − (cid:19) , (12)where λ ∗ a , can be obtained by inversion of (11), once theglobal arrival rate of molecules is known.It is worth to note that although the system model illustratedin Fig. 1 assumes that the target of the drug is modeled as aspherical bio-nanomachine, in order to resemble some typesof cells, the application of the proposed model is much moregeneral. In more detail, the target of the drug can be modeledas a surface with a generic shape provided with receptors.If fact, as it appears from (11) and (12), the per-receptorassimilation and rejection rates are functions of the receptordistance from the center of transmission nanomachine, and canbe easily evaluated once a measure λ ∗ a is available. If we assume that receptors are uniformly distributed onthe RX surface, and by using the symmetry around the axisconnecting the center of TX with the center of RX, it ispossible to find more results. By partitioning the RX surfaceinto R RX areas, it results that the surface surrounding eachreceptor has an average size S = π r RX / R RX . The angle α = β corresponding to the receptor located at distance d + r RX willbe equal to β = / √ R RX due to geometrical considerations. Ifwe divide the remaining surface of half semi-sphere (that forpositive value of α , the same reasoning holds for negativeones) in a number of spherical zones with angle ∆ α , allreceptors in these zones have the same distance from thecenter of TX, thus they see the same local concentration ofdrug molecules, and thus the same arrival rate λ ∗ a , j . If weassume that the area S dedicated to each receptor has roughly asquare shape, it results that S = ( r RX ∆ α ) , thus ∆ α = √ S / r RX = √ π/ R RX , and the number of spherical zones F are given by F = π − β ∆ α . Since the area of each k -th spherical zone identifiedby α k = k ∆ α is given by S k = π r RX ( sin ( k ∆ α ) − sin (( k − ∆ α ))(see also Fig. 1), the number of receptors in such a zone is n k = S k S = R RX (cid:16) sin ( k √ π/ R RX ) − sin (( k − √ π/ R RX ) (cid:17) λ ∗ r = λ ∗ a , ( d + r RX ) F + X k = − F − n k r ( α k ) (cid:18) r ( α k ) λ ∗ a , ( d + r RX ) T traf f − (cid:19) , (14)and λ ∗ a , can be obtained by λ ∗ a , = λ ∗ a ( d + r RX ) P F + k = − F − n k r ( α k ) . (15)In order to consider a di ff erent distribution of receptors, byassuming that the symmetry with respect to d still holds, it isenough to change the values of n k , given that the condition P F + k = − F − n k = R RX is satisfied. C. Application to drug delivery
In this sub-section we illustrate the application of ourcongestion model to drug delivery systems. Without loss ofgenerality, we refer to the occupancy theory [59], since it isthe most commonly accepted theory. It says that the magnitudeof drug response depends on the occupancy ratio of receptorsby drug molecules. What we illustrate in what follows cansimilarly be adapted to the rate theory.As shown in [61], in some cases (e.g. competitive antagonistdrugs) a full drug response can be produced even at lowreceptor occupancy, and it results that excess receptors, notbound to drug molecules, are no further needed to obtainthe maximum drug response. Thus, in this case study weassume the usage of a drug type which produces the desiredresponse when at least a fraction f of receptors of the targetcells are bound to the drug molecules. Thus, in order tosuitably configure the drug delivery system, the first stepconsists of determining the minimum drug release rate thatallows achieving the desired percentage of bound receptors.By using the results shown in [15], it is easy to find that a good approximation of the arrival rate value λ o = λ ∗ a + λ ∗ r isgiven by λ o = λ ∗ a + λ ∗ r = r RX d R RX r r , rx π r RX + R RX r r , rx Q ∆ t . (16)Thus, by applying (14) and (15) in (16), it is easy to finda relationship between the release rate and λ ∗ a . In addition, byusing the Little’s law [64], it is immediate to find that the meannumber of busy receptors (servers) in the system is given by λ ∗ a T tra f f , and the desired fraction f results equal to the overallutilization coe ffi cient ρ , which is equal to ρ = λ ∗ a T tra f f R RX = f . (17)IV. P erformance E valuation The performance of the system has been evaluated by usingthe BiNS2 simulator [66], [67], which is a Java packagedesigned to simulate nano-scale biological communications in3D. The approach of BiNS2 is fine grained: the position ofeach element, nanomachine or molecule, is evaluated at eachsimulation step, and collisions are managed according to apartially inelastic model, described by means of the coe ffi cientof restitution [68]. Not only molecules, but also nanomachinescan be either fixed or mobile. The surrounding environmentcan be unbounded or bounded. In the latter case, the boundingsurface can have di ff erent shapes; currently the supported onesare sphere, cylinder, and cube, or a combination of them [51].In addition, BiNS2 allows modeling a receiver node with afinite number of receptors R RX . Each receptor implements afinite, non-negligible reception time (tra ffi cking time). When amolecule hits one of the molecule-compliant receptors, if it isnot busy in another bond, it locks the relevant receptor for anamount of time whose distribution of the can be selected froma number of known statistics. In this work, the tra ffi cking timeis exponentially distributed with average value T tra f f , whereasthe transmission rate at TX is fixed for the whole simulationduration. Instead, if the receptor is busy, the molecule isbounced back as it would have hit a portion of the surfaceof the receiver without receptors. In the scenario of this paper,nanomachines are fixed, whereas molecules move accordingto the Brownian motion and as results of collisions, whichcan occur with nanomachines or among themselves.The simulation reliability of BiNS2 has been assessed bytuning simulation results with lab experiments [67]. Finally,the BiNS2 package implements an octree-based computationapproach [67], which uses a dynamic splitting of the simulatedenvironment into cubes of di ff erent size in order to parallelizethe simulation, so as to benefit of the multi-thread capabilitiesof modern multi-core computer architectures, and thus stronglyreducing the simulation time. In the simulation of the scenarioconsidered in this paper, the simulation environment is un-bounded; however, we implemented a cubic virtual boundarywith side 1 mm to implement the octree algorithm. TheTX nanomachine is located at the center of the cube. Whenmolecules exit the cube, they are removed from the simulation,in order to limit the computational burden of the simulation.From the above description, it appears that BiNS2 implements TABLE I: Simulation parameters Symbol Description Value T Temperature 310 K e Coe ffi cient of restitution (part. inelastic collisions) 0.95 η Viscosity 0.0011 Pl [68] r RX Radius node RX 2.5 µ mr TX Radius node TX 2.5 µ mR RX Amount of surface receptors (node RX) 10000 r c , tx Radius emitted molecules 1.75 nm r r , rx Receptor radius (RX) 4 nm T traf f Tra ffi cking time 2 or 4 s ∆ t Emission period (TX) 20 ms d RX-TX distance 26.5 µ m
40 50 60 70 80 90 100405060708090100 Simulation data S y n t he t i c da t a Fig. 2: Q-Q plot of arrivals (both absorbed and rejectedmolecules in intervals of 0.1 s) obtained via simulation versussynthetic Poisson data with the same average.exactly the receiver model that we consider in this paper. Themain simulation parameters, together with their descriptionsand values, are reported in Table I. In the next subsection,we present the simulation results of the scenario described insection III.
A. Numerical results
Before evaluating system performance, it is necessary toverify the suitability of the assumption of Poisson distributedarrivals. Fig. 2 presents the Q-Q plot of simulation dataversus synthetic Poisson data generated with the same averagevalue. The simulation data include both absorbed and rejectedmolecules, gathered in intervals of 0.1 seconds. The resultingagreement is evident.From the results in sections III, it can be argued that theconcentration of molecules is higher on the portion of thereceiver surface facing the TX. In this portion of the surface,there are a lot of molecules which cannot establish a bondwith receptors since most of them are already busy. Instead,the concentration of molecules close to the opposite portionof the receiver surface is much lower, since the molecules,following the negative gradient of concentration, tend to leavethe RX node. As shown in section III-B, this phenomenonhas a significant impact also on assimilations (and rejections).This observation is confirmed in Fig. 3, which shows the mapof assimilations for a tra ffi cking time equal to T ta f f = ff erent times, corresponding to di ff erent amount of assimi-lations. N a represents the cumulative number of assimilations for each receptor, and it is plotted by means of a discretevalue, chromatic scale for all the R RX receptors of RX. Eachreceptor is identified by its spherical coordinates, that is theazimuth θ , whose range is (- π , π ] radians, and the altitude φ ,whose range is (- π/ π/
2] radians. The direction connectingthe center of TX and RX corresponds to ( θ, φ ) = (0 , N a values),in all considered samples, is located in the region around thecenter of the plot, that is ( θ, φ ) = (0 , θ and φ close to zero continue remainingclose to the surface of the RX, due to negative gradient of theconcentration, which “pushes” them towards large values of d . This facilitates the creation of ligand-receptor bonds on alarger area. However, in the region with coordinates close to( θ, φ ) = (0 , N a are often more than twice thevalue in the remaining of the RX surface. Since these receptorsof the RX receive a larger portion of binding attempts, theyare also responsible for most of rejections. We now showsimulation results obtained by using both the proposed modeland the one illustrated in [16], suitably adapted to calculatethe rejection rate, as described in section III-A. In order tocompare simulations with the theoretical models, we haveassumed that they have the same mean arrival rate λ a = λ ∗ a (and thus the same mean number of busy receptors). We haveevaluated the rejection rate, since λ a is the only parameterthat the RX can estimate reliably. In the simulation, we havetracked the number of rejection events, whereas the theoreticalvalue in (9) has been determined by evaluating the stateprobabilities numerically. As for the proposed model, the rateof rejection (i.e. unabsorbed drug molecules due to collisionswith busy receptors) is evaluated by means of (14). In addition,in the simulation we have verified that the ratio between thecollisions of molecules with any of the deployed receptors,either busy or not, and those with the portion of surface of theRX not covered by receptors, closely match the theoreticalvalue of R RX ( r r , rx + r c , rx ) r RX , where r r , rx is the receptor radius, r c , tx is the molecules radius, and r RX is the RX node radius.Results are shown in Fig. 4. In Fig. 4.a, the ordinate axisreports the rejection rate for the case T tra f f = T tra f f = λ r , which means that it is not able to e ff ectivelydetect congestion. Instead, our proposed model exhibits anexcellent match with simulation results in all cases. Thisis an expected result, since in [16] it is assumed that theconcentration of molecules is the same all over the surface ofthe receiver. Instead, it could be not true, especially for lowtransmission rates and for short tra ffi cking times. In addition, and this is the most important reason, receptors do not behaveas a pool of servers, which can be invoked and used upon anew user enters the system (M / M / m / m model). Instead, theyare isolated service elements, and this nature must be takeninto account. Even if the symmetric model tries to mitigatethis e ff ect, by using a state dependent arrival process, it is notsu ffi cient. In particular, at low T tra f f values, the di ff erencebetween the results of the symmetric model and the simu-lation increases, since it is expected that the RX is eager toassimilate drug molecules and to re-present available receptors.Nevertheless, this is true only in unlikely case of uniformconcentration. Instead, an interesting behavior appears in thecase T tra f f =
4s when the absorption rate is very high, largerthan 1000 molecules / s. In this case, as shown in Fig. 4.a, theestimation of the rejection rate given by the symmetric modeltends to converge to the simulation values, which is also thevalue estimated by our model. Also this phenomenon can beeasily explained. In fact, when the T tra f f value increases, thereceiver is less prompt to free receptors. This means that whenthe absorption rate increases beyond a threshold value, due toa significant increase of the concentration of drug moleculesnearby RX, even if the concentration close to the receiversurface is not homogeneous, it is so high that the limitingfactor becomes the number of free receptors (about the 50%in average), and thus the di ff erences between the simulationand the two theoretical models tends to vanish. This trend isonly barely visible for the case of T tra f f = s , due to the factthat, at the maximum transmission rate, the average numberof busy receptors is still low (about 32%). In all the othercases, the di ff erence between the symmetric and the proposedmodel is in the range of at least one order of magnitude.Fig. 5 shows the blocking probability, evaluated as λ ∗ r λ o , where λ o = λ ∗ r + λ ∗ a is the system arrival rate. The abscissa reports theabsorption rate, expressed in molecules per second. The nete ff ect of the behaviors observed in the previous Fig. 4 is anearly liner increase of the blocking probability, as a functionof the arrival rate, which is typical of pure loss systems inoverload. Finally, Fig. 6 shows the e ff ect of varying the totalnumber of receptors R RX by 20% on the rejection rate, whilekeeping constant λ ∗ a , indicated in the abscissa. As expected,when R RX increases, the rejection rate also decreases. Instead,when R RX decreases, a marked increase of the rejection rateappears. Again, this phenomenon can be explained with thefact that receptors are not a multiplexed resource, but they arerather multiple, single resources. By increasing the distancebetween TX and RX, we have also investigated the occurrenceof any marked decrease in the rejection rate, due to the factthat di ff erent concentration values should be less significantfor large d (see also Fig. 3). However, we have found a nearlynegligible decrease for the rate of rejected molecules with d for T tra f f = s , not reported in this paper. This is probably dueto the fact that the phenomenon visible in Fig. 3 is due to notonly to the transmission distance, but also to the disturbinge ff ect of the presence of the node RX itself. This is anotherconfirmation that, when the system absorption rate λ ∗ a is keptconstant, the factor with the most significant impact is thenumber of receptors, since it is not a multiplexed resource. −4 −3 −2 −1 0 1 2 3 4−2−1.5−1−0.500.511.52 A l t i t ude φ (r ad ) Azimut θ (rad)a) Total number of assimilations: 7850 −4 −3 −2 −1 0 1 2 3 4−2−1.5−1−0.500.511.52 A l t i t ude φ (r ad ) Azimut θ (rad)b) Total number of assimilations: 48671−4 −3 −2 −1 0 1 2 3 4−2−1.5−1−0.500.511.52 A l t i t ude φ (r ad ) Azimut θ (rad)c) Total number of assimilations: 72891 −4 −3 −2 −1 0 1 2 3 4−2−1.5−1−0.500.511.52 A l t i t ude φ (r ad ) Azimut θ (rad)d) Total number of assimilations: 117920 1 ≤ N a <5 5 ≤ N a <10 10 ≤ N a <15 15 ≤ N a <20 N a ≥ Fig. 3: Maps of the number of assimilations for each surface receptor ( N a ), as a function of its spherical coordinates, azimuthand the altitude, for increasing numbers of total assimilation: a) 7850, b) 48761, c) 72891, and d) 117920. B. Lesson learned towards drug delivery
In the previous subsection, we have analyzed the proposedmodel through simulations. Now, we highlight some keyissues learned from this analysis, and focus on the operationalprocedure for designing an e ff ective, localized drug deliverysystem.First, in order to make use of the proposed model, it isnecessary to characterize the underlying system. This can bedone by means of lab experiments, or by using results alreadypresent in the literature. The goal of these experiments consistsof estimating a small number of parameters, necessary forcharacterizing the system behavior. These parameters are theaverage tra ffi cking time T tra f f , and the number of surfacereceptors of the considered type, R RX . Once the values ofthese two parameters are known, and the TX-RX distance d is known, it is possible to make use of the model for anyvalue of the emission rate Q ∆ t in order to evaluate the averagenumber of busy receptors. Finally, according to the drug used,it is necessary to know the minimum fraction of receptors to bebound to drug molecules in order to maximize the drug e ff ects. As mentioned above, for some types of drug this fraction couldbe quite low. For instance, in [61], the author shows that oftenonly 5-10% occupancy is needed to produce a full responsewhen agonist drugs are used. For each drug, this value can bederived by in-vitro experiments. Thus, the final step consist ofusing the results obtained in subsection III-C for estimatingthe optimal drug release rate Q ∆ t , which guarantees the desiredfraction of receptors bound to drug molecules, that is at least f , without drug overloading. Fig. 7.a shows the percentageof busy receptors as a function of the release rate of drugmolecules, for the tra ffi cking time values already used in thispaper, 2s and 4s, respectively. As expected, the value of f increases with the rate Q ∆ t . This sublinear increase is due tocongestion. Fig. 7.b depicts the minimum release rate of drugmolecules necessary to achieve the target occupancy f , versusthe transmission range, for di ff erent values of the tra ffi ckingtime. A first comment is that the dependency of drug releaserate is nearly linear with the transmission range. In addition,the tra ffi cking time has a significant influence on the slope ofthe release rate curve. In particular, the lower the tra ffi ckingtime, the prompter is the target cell to absorb the drug, which −2 Absorption rate (molecules/s) R e j e c t i on r a t e ( m o l e c u l e s / s ) a) T traff =4 s SimulationsProposed modelSymmetric model10 −2 Absorption rate (molecules/s) R e j e c t i on r a t e ( m o l e c u l e s / s ) b) T traff =2 s SimulationsProposed modelSymmetric model Fig. 4: Rejections as a function of the absorption rate: a) T tra f f equal to 4 seconds, and b) T tra f f equal to 2 seconds. B l o ck i ng p r ob . T traff =1/ µ =2 sT traff =1/ µ =4 s Fig. 5: Blocking probability as a function of the absorptionrate.makes it necessary to transmit more molecules to maintain thetarget occupancy f . V. C onclusion In this paper, we have introduced the concept of congestionin di ff usion-based molecular communications, and proposeda model for both illustrating the dynamical behavior ofthe phenomenon and analyzing the root causes of it. Theconsidered continuous emission of molecules is typical oflocalized drug delivery systems. Our proposed model of thereceiver nanomachine includes an M / M / /
500 550 600 650 700 750 800100200300400500600 Absorption rate (mol/s) R e j e c t i on r a t e ( m o l / s ) a) T traff =4 s R RX =10000R RX =12000R RX =8000500 550 600 650 700 750 800050100150200250 Absorption rate (mol/s) R e j e c t i on r a t e ( m o l / s ) b) T traff =2 s R RX =10000R RX =12000R RX =8000 Fig. 6: Rejection rate as a function of the absorption rate fordi ff erent number of receptors and tra ffi cking times. F r a c t i on o f bu sy r e c ep t o r s a) T traff =2sT traff =4s10 20 30 40 50 60010203040 TX−RX distance (um) E m i ss i on r a t e ( t hou s and s o f m o l e c u l e s / s ) b) T traff =2s, f=0.1T traff =2s, f=0.2T traff =4s, f=0.1T traff =4s, f=0.2 Fig. 7: Percentage of busy receptors as a function of drugmolecules release rate.considered a set of multiplexed resources such as a server pool.They behave as isolated, multiple nano-receivers, independentof each other. An extensive numerical analysis shows thatthis model is e ff ective in modeling congestion conditions,especially in the range of distances typical of di ff usion-basedmolecular communication.The results of this work can be used for implementing ratecontrol algorithms for molecular communications, in particular for drug delivery systems. We have detailed the operationalprocedure for determining the suitable drug delivery rate froma set of implanted emitting nanomachines, which release drugnanoparticles close to the target cells.In addition, since the arrival rate can be measured by thereceiver bio-nanomachine by counting the absorbed molecules,the system can also be designed with adaptive features. Inparticular, it is possible to add a further controller nanoma-chine, which can estimate whether the used transmission rateis optimal, and adapt its value to any change of the operationalconditions. Release rate adjustments can be triggered by usingfeedback messages [8] sent when the observed deviations fromthe original design conditions are significant.The drug delivery system analyzed in this paper is analyzedin a quasi-static environment, such as the extracellular matrix.Drug delivery in more dynamic environments, such as bloodvessels, we will be the object of future investigations. In suchscenarios di ff erent configurations can be studied, such as fixedtransmitters, anchored to vessel wall, and mobile targets (e.g.circulating tumor cells), or mobile transmitters and mobiletargets, or mobile transmitters and fixed targets (e.g. solidtumors). R eferences [1] I. F. Akyildiz, F. Brunetti, and C. Blázquez, “Nanonetworks: A newcommunication paradigm,” Comput. Netw. , vol. 52, no. 12, pp. 2260–2279, Aug. 2008.[2] T. Nakano, M. Moore, F. Wei, A. Vasilakos, and J. Shuai, “Molecularcommunication and networking: Opportunities and challenges,”
IEEETransactions on NanoBioscience , vol. 11, no. 2, pp. 135–148, 2012.[3] T. M. Allen and P. R. Cullis, “Drug delivery systems: Entering themainstream,”
Science , vol. 303, no. 5665, pp. 1818 – 1822, 2004.[4] K. Winstein and H. Balakrishnan, “TCP ex machina: Computer-generated congestion control,” in
ACM SIGCOMM’13 , 2013, pp. 123–134.[5] A. Sivaraman, K. Winstein, Thaker, and H. Balakrishnan, “An ex-perimental study of the learnability of congestion control,” in
ACMSIGCOMM’14 , 2014, pp. 479–490.[6] V. Jacobson, “Congestion avoidance and control,” in
ACM SIGCOMM’88 . New York, NY, USA: ACM, 1988, pp. 314–329.[7] K. K. Ramakrishnan and R. Jain, “A binary feedback scheme forcongestion avoidance in computer networks,”
ACM Trans. Comput. Syst. ,vol. 8, no. 2, pp. 158–181, May 1990.[8] T. Nakano, Y. Okaie, and A. V. Vasilakos, “Transmission rate controlfor molecular communication among biological nanomachines,”
IEEEJournal on Selected Areas in Communications , vol. 31, no. 12, suppl.,2013.[9] L. Felicetti, M. Femminella, G. Reali, T. Nakano, and A. V. Vasilakos,“TCP-like molecular communications,”
IEEE Journal on Selected Areasin Communications , Dec. 2014.[10] J. Philibert, “One and a half century of di ff usion: Fick, einstein, beforeand beyond,” Di ff usion Fundamamentals , vol. 4, pp. 6.1–6.19, 2006.[11] D. Lau ff enburger and J. Linderman, Receptors: Models for Binding,Tra ffi cking, and Signalling, . Oxford University Press, 1996.[12] P. Barta, H. Bjorkelund, and K. Andersson, “Circumventing the require-ment of binding saturation for receptor quantification using interactionkinetic extrapolation,” Nucl Med Commun , vol. 32, no. 9, pp. 863–867,Sep 2011.[13] D. J. Freeman, K. McDorman, S. Ogbagabriel, C. Kozlosky, B. B. Yang,S. Doshi, J. J. Perez-Ruxio, W. Fanslow, C. Starnes, and R. Radinsky,“Tumor penetration and epidermal growth factor receptor saturation bypanitumumab correlate with antitumor activity in a preclinical model ofhuman cancer,”
Mol. Cancer , vol. 11, p. 47, 2012.[14] H. Yilmaz, A. Heren, T. Tugcu, and C.-B. Chae, “Three-dimensionalchannel characteristics for molecular communications with an absorbingreceiver,”
Communications Letters, IEEE , vol. 18, no. 6, pp. 929–932,June 2014. [15] A. Akkaya, H. Yilmaz, C. Chae, and T. Tugcu, “E ff ect of receptordensity and size on signal reception in molecular communication viadi ff usion with an absorbing receiver,” Communications Letters, IEEE ,vol. 19, no. 2, pp. 155–158, Feb 2015.[16] M. Pierobon and I. Akyildiz, “Noise analysis in ligand-binding receptionfor molecular communication in nanonetworks,”
IEEE Transactions onSignal Processing , vol. 59, no. 9, pp. 4168 –4182, September 2011.[17] B. Koo, H. Yilmaz, C.-B. Chae, and E. A., “Detection algorithms formolecular MIMO,” in
IEEE ICC’15 , 2015.[18] M. F. L. Felicetti and G. Reali, “Smart antennas for di ff usion-basedmolecular communications,” in ACM Nanocom 2015 , 2015.[19] M. D.A., “Stochastic approach to chemical kinetics,”
J Appl Probab ,vol. 4, pp. 413–478, 1967.[20] T. Nakano, T. Suda, Y. Okaie, M. Moore, and A. Vasilakos, “Molecularcommunication among biological nanomachines: A layered architectureand research issues,”
NanoBioscience, IEEE Transactions on , vol. 13,no. 3, pp. 169–197, Sept 2014.[21] A. Aijaz, A. Aghvami, and M. Nakhai, “On error performance ofnetwork coding in di ff usion-based molecular nanonetworks,” Nanotech-nology, IEEE Transactions on , vol. 13, no. 5, pp. 871–874, Sept 2014.[22] V. Loscri, C. Marchal, N. Mitton, G. Fortino, and A. Vasilakos, “Securityand privacy in molecular communication and networking: Opportunitiesand challenges,”
NanoBioscience, IEEE Transactions on , vol. 13, no. 3,pp. 198–207, Sept 2014.[23] M. Moore, T. Suda, and K. Oiwa, “Molecular communication: Modelingnoise e ff ects on information rate,” IEEE Transactions on NanoBio-science , vol. 8, no. 2, pp. 169–180, 2009.[24] T. Nakano, Y. Okaie, and J.-Q. Liu, “Channel model and capacityanalysis of molecular communication with brownian motion,”
Commu-nications Letters, IEEE , vol. 16, no. 6, pp. 797–800, 2012.[25] M. . Kuran, H. B. Yilmaz, T. Tugcu, and B. Ozerman, “Energy model forcommunication via di ff usion in nanonetworks,” Nano CommunicationNetworks , vol. 1, no. 2, pp. 86 – 95, 2010.[26] S. Kadloor, R. Adve, and A. Eckford, “Molecular communication usingbrownian motion with drift,”
IEEE Transactions on NanoBioscience ,vol. 11, no. 2, pp. 89–99, June 2012.[27] K. V. Srinivas, A. Eckford, and R. Adve, “Molecular communicationin fluid media: The additive inverse gaussian noise channel,”
IEEETransactions on Information Theory , vol. 58, no. 7, pp. 4678–4692,2012.[28] B. Atakan and O. B. Akan, “On molecular multiple-access, broadcast,and relay channels in nanonetworks,” in
ICST BIONETICS 2010 , 2010.[29] ——, “Deterministic capacity of information flow in molecular nanonet-works,”
Nano Communication Networks , vol. 1, no. 1, pp. 31 – 42, 2010.[30] M. U. Mahfuz, D. Makrakis, and H. T. Mouftah, “On the characterizationof binary concentration-encoded molecular communication in nanonet-works,”
Nano Communication Networks , vol. 1, no. 4, pp. 289 – 300,2010.[31] M. Pierobon and I. Akyildiz, “A physical end-to-end model for molec-ular communication in nanonetworks,”
IEEE Journal on Selected Areasin Communications , vol. 28, no. 4, pp. 602–611, 2010.[32] ——, “Di ff usion-based noise analysis for molecular communication innanonetworks,” IEEE Transactions on Signal Processing , vol. 59, no. 6,pp. 2532–2547, 2011.[33] T. Nakano and J.-Q. Liu, “Design and analysis of molecular relaychannels: An information theoretic approach,”
IEEE Transactions onNanoBioscience , vol. 9, no. 3, pp. 213–221, 2010.[34] T. Nakano and J. Shuai, “Repeater design and modeling for molecularcommunication networks,” in
IEEE INFOCOM Workshop , 2011, pp.501–506.[35] N. Farsad, A. Eckford, and S. Hiyama, “A markov chain channel modelfor active transport molecular communication,”
Signal Processing, IEEETransactions on , vol. 62, no. 9, pp. 2424–2436, May 2014.[36] A. Guney, B. Atakan, and O. Akan, “Mobile ad hoc nanonetworkswith collision-based molecular communication,”
IEEE Transactions onMobile Computing , vol. 11, no. 3, pp. 353–366, March 2012.[37] H. ShahMohammadian, G. G. Messier, and S. Magierowski, “Optimumreceiver for molecule shift keying modulation in di ff usion-based molec-ular communication channels,” Nano Communication Networks , vol. 3,pp. 183–195, 2012.[38] B. Atakan and O. Akan, “Deterministic capacity of information flow inmolecular nanonetworks,”
Nano Communication Networks , vol. 1, pp.31–42, 2010.[39] W. Guo, S. Wang, A. Eckford, and J. Wu, “Reliable communicationenvelopes of molecular di ff usion channels,” Electronics Letters , vol. 49,no. 19, pp. 1248–1249, 2013. [40] I. Llatser, A. Cabellos-Aparicio, M. Pierobon, and E. Alarcon, “De-tection techniques for di ff usion-based molecular communication,” IEEEJournal on Selected Areas in Communications , vol. 31, no. 12, supple-ment, 2013.[41] L.-S. Meng, P.-C. Yeh, K.-C. Chen, and I. Akyildiz, “On receiver designfor di ff usion-based molecular communication,” Signal Processing, IEEETransactions on , vol. in press, 2014.[42] M. Mahfuz, D. Makrakis, and H. Mouftah, “A comprehensive studyof sampling-based optimum signal detection in concentration-encodedmolecular communication,”
NanoBioscience, IEEE Transactions on ,vol. 13, no. 3, pp. 208–222, Sept 2014.[43] M. Pierobon and I. Akyildiz, “A statistical-physical model of inter-ference in di ff usion-based molecular nanonetworks,” Communications,IEEE Transactions on , vol. 62, no. 6, pp. 2085–2095, June 2014.[44] D. Kilinc and O. Akan, “Receiver design for molecular communication,”
Selected Areas in Communications, IEEE Journal on , vol. 31, no. 12,pp. 705–714, December 2013.[45] H. B. Yilmaz and C.-B. Chae, “Simulation study of molecular com-munication systems with an absorbing receiver: Modulation and ISImitigation techniques,”
Simulation Modelling Practice and Theory ,vol. 49, pp. 136 – 150, 2014.[46] S. EL Andaloussi, I. Mager, X. O. Breakefield, and M. J. Wood,“Extracellular vesicles: biology and emerging therapeutic opportunities,”
Nat Rev Drug Discov , vol. 12, no. 5, pp. 347–357, May 2013.[47] D. C. Gadsby, “Ion channels versus ion pumps: the principal di ff erence,in principle,” Nat. Rev. Mol. Cell Biol. , vol. 10, no. 5, pp. 344–352, May2009.[48] L. Felicetti, M. Femminella, G. Reali, J. Daigle, M. Malvestiti, andP. Gresele, “Modeling CD40-based molecular communications in bloodvessels,”
IEEE Transactions on NanoBioscience , vol. 13, no. 3, pp. 230–243, 2014.[49] D. R. Elias, A. Poloukhtine, V. Popik, and A. Tsourkas, “E ff ect ofligand density, receptor density, and nanoparticle size on cell targeting,” Nanomedicine: Nanotechnology, Biology and Medicine , vol. 9, no. 2,pp. 194 – 201, 2013.[50] J. Tan, A. Thomas, and Y. Liu, “Influence of red blood cells onnanoparticle targeted delivery in microcirculation,”
Soft Matter , vol. 8,pp. 1934–1946, 2012.[51] L. Felicetti, M. Femminella, and G. Reali, “Simulation of molecularsignaling in blood vessels: Software design and application to athero-genesis,”
Nano Communication Networks , vol. 4, no. 3, pp. 98 – 119,2013.[52] L. L. M. Heaton, E. López, P. K. Maini, M. D. Fricker, and N. S.Jones, “Advection, di ff usion, and delivery over a network,” Phys. Rev.E , vol. 86, p. 021905, Aug 2012.[53] Y. Chahibi, M. Pierobon, S. O. Song, and I. Akyildiz, “A molecularcommunication system model for particulate drug delivery systems,”
Biomedical Engineering, IEEE Transactions on , vol. 60, no. 12, pp.3468–3483, Dec 2013.[54] Y. Chahibi, M. Pierobon, and I. Akyildiz, “Pharmacokinetic modelingand biodistribution estimation through the molecular communicationparadigm,”
Biomedical Engineering, IEEE Transactions on , vol. PP,no. 99, pp. 1–1, 2015.[55] Y. Chahibi, I. Akyildiz, S. Balasubramaniam, and Y. Koucheryavy,“Molecular communication modeling of antibody-mediated drug deliv-ery systems,”
Biomedical Engineering, IEEE Transactions on , vol. 62,no. 7, pp. 1683–1695, July 2015.[56] D. Colquhoun, “The quantitative analysis of drug-receptor interactions:a short history,”
Trends in Pharmacological Sciences , vol. 27, no. 3, pp.149 – 157, 2006, 75 years of the British Pharmacological Society.[57] T. Kenakin, “Principles: Receptor theory in pharmacology,”
Trends inPharmacological Sciences , vol. 25, no. 4, pp. 186 – 192, 2004.[58] H. P. Rang, “The receptor concept: pharmacology’s big idea,”
Br. J.Pharmacol. , vol. 147 Suppl 1, pp. 9–16, Jan 2006.[59] A. Clark,
The Mode of Action of Drugs on Cells . Eduard Arnold andCo., 1933.[60] W. D. M. Paton, “A theory of drug action based on the rate of drug-receptor combination,”
Proceedings of the Royal Society of London B:Biological Sciences , vol. 154, no. 954, pp. 21–69, 1961.[61] D. Lambert, “Drugs and receptors,”
Continuing Education in Anaesthe-sia, Critical Care and Pain , vol. 4, no. 6, pp. 181–184, 2004.[62] D. R. O’Brien, “The inverse relationship between drug a ffi nity ande ff ectiveness: Prediction under rate theory, paradox under occupancytheory,” Medical Hypotheses , vol. 23, no. 3, pp. 327 – 333, 1987.[63] A. Einolghozati, M. Sardari, and F. Fekri, “Design and analysis of wire-less communication systems using di ff usion-based molecular communi- cation among bacteria,” Wireless Communications, IEEE Transactionson , vol. 12, no. 12, pp. 6096–6105, December 2013.[64] L. Kleinrock,
Theory, Volume 1, Queueing Systems . Wiley-Interscience,1975.[65] S. Kuo, D. Hammer, and D. Lau ff enburger, “Simulation of detachmentof specifically bound particles from surfaces by shear flow,” BiophysicalJournal , vol. 73, no. 1, pp. 517 – 531, 1997.[66] L. Felicetti, M. Femminella, and G. Reali, “A simulation tool fornanoscale biological networks,”
Nano Communication Networks , vol. 3,no. 1, pp. 2–18, 2012.[67] L. Felicetti, M. Femminella, G. Reali, P. Gresele, and M. Malvestiti,“Simulating an in vitro experiment on nanoscale communications byusing bins2,”
Nano Communication Networks , vol. 4, no. 4, pp. 172 –180, 2013.[68] D. Gidaspow and J. Huang, “Kinetic theory based model for blood flowand its viscosity,”