# Triangulation sensing

Triangulation sensing is a theory describing the computational steps of a cell containing on its surface small windows, to recover the location of a source emitting random particles in a medium. This is particularly relevant for neuron navigation in the Brain [1][2].Reconstructing the source location allows the cell to triangulate its position. The reconstruction steps of the gradient source from the fluxes of diffusing particles arriving to small absorbing receptors are 1-arrival of the Brownian particles to the small windows 2- counting of the particles 3-inversion of the Laplace's equation to estimate the position from the fluxes 4- possible noise reduction by applying the same procedure to several triplets [3] .

## Computing the diffusion fluxes to small window with respcet to the source location

The mathematical formulation consists inconsidering diffusing cues that have to bind to N narrow windows located on the surface of a three dimensional shaped object, typically a ball (in dimension 3) or a disk in dimension 2. For a ball B(a) of radius a. M can range between 10 and 50, but can be much higher for other receptor types. Individual cue molecules are described as Brownian particles. The cues are released from a source at position $x_0$ outside the ball. The triangulation sensing method consists in reconstructing the source position from estimated steady-state fluxes at each narrow window for fast binding (i.e. the probability density has an absorbing boundary condition at the windows).[4]

The steady-state gradient $P_0(x)$ is obtained by integrating the diffusion equation from 0 to infinity, which is equivalent to resetting a particle after it disappears through a window. It is given by the solution of the mixed-boundary value problem (where D is the diffusion coefficient)

$D\Delta P_0(x) = -Q \delta(x-x_0) \hbox{ for } x\, \in \,\R^3 -\Omega$

$\frac{\partial P_0}{\partial n}(x) = 0\hbox{ for } x, \in\, \partial\Omega_r$

$P_0(x) = 0 \hbox{ for } x \in \partial \Omega_a.$

The probability fluxes associated with $P_0(x)$ at each individual window $\Omega_a(\epsilon_j)$ can be computed from the fluxes and depend on the specific window arrangement and the domain $\Omega$. The parameter Q>0 can be calibrated so that there are a fixed number of particles located in a given volume. At infinity, the density $P_0(x)$ has to tend to zero in three dimensions. More complex domains could be studied if their associated Green's function can be found explicitly.

Note that in dimension, although the probability density $P_0(x)$ diverges when $x\rightarrow\infty$, the splitting probability between windows is finite because it is the ratio of the steady-state flux at each window divided by the total flux through all windows:

$J_k= \frac{ \int_{\Omega_{k}} \frac{\partial P_0(x)}{\partial n} dS_{x}}{\sum_{q} \int_{\partial \Omega_{q}} \frac{\partial P_0(x)}{\partial n} dS_{x}}.$

In two-dimensions, due to the recurrent property of the Brownian motion, the probability to hit a window before going to infinity is one, thus the total flux is one:

$\sum_{q} \int_{\partial\Omega_{q}} \frac{\partial P_0(x)}{\partial n} dS_{x}=1.$

The construction of the solution uses an inner and outer solution. The inner solution is constructed near each small window [5] by scaling the Arc length s and the distance to the boundary $\eta$ by $\bar \eta=\frac{\eta}{\epsilon}$ and $\bar s=\frac{s}{\epsilon}$, so that the inner problem reduces to the classical two-dimensional Laplace equation[6]

$\Delta w=0 \hbox{ in } \R_+^{2}, \frac{\partial w}{\partial n}=0 \hbox{ for } |\bar s|>\frac{1}{2}, \bar \eta=0 , w(\bar s, \bar \eta)=0 \hbox{ for } |\bar s|<\frac{1}{2}, \bar \eta=0.$

The far field behavior for $|x|\rightarrow \infty$ and for each hole i=1, 2 is $w_i(x) \approx A_i \{\log|x-x_i|-\log \epsilon\} +o(1),$ where A_i is the flux $A_i=\frac{2}{\pi} \int_0^{1/2} \frac{\partial w(0,\bar s)}{\partial \bar \eta}d\bar s.$

The general solution of equation for n=2 is obtained from the outer solution of the external Neumann-Green's function

$-\Delta_{x} G(x,y) = \delta(x-y) \hbox{ for } x \in \,\R_{+}^2$

$\frac{\partial G}{\partial n_{x}}(x,x_0) = 0\hbox{ for } x \in {\partial\R_{+}^2}.$

Given for $x, x_0 \in \R_{+}^2$ by $G(x,x_0)=\frac{-1}{2\pi }\left(\ln |x-x_0| + \ln \left|x-\bar{x_0}\right| \right),$ where $\bar{x_0}$is the symmetric image of $x_0$through the boundary axis 0z. The uniform solution is the sum of inner and outer solution (Neuman-Green's function) $P(x,x_0)=G(x,x_0) +A_1 \{\log|x-x_1|-\log \epsilon\}+A_2 \{\log|x-x_2|-\log \epsilon\} +C,$ where $A_1,A_2,C$ are constants to be determined. To that purpose, we study the behavior of the solution near each point $x_i$. In the boundary layer, we get

$P(x,y) \approx A_i \{\log|x-x_i|-\log \epsilon\}.$

Using this condition on each window, we obtain the two conditions:

$G(x_1,x_0)+A_2 \{\log|x_1-x_2|-\log \epsilon\} +C=0$

$G(x_2,x_0) +A_1 \{\log|x_2-x_1|-\log \epsilon\} +C=0.$

Due to the recursion property of the Brownian motion in dimension 2, there are no fluxes at infinity, thus the conservation of flux gives: $\int_{\partial\Omega_{1}} \frac{\partial P(x,y)}{\partial n}dS_{x} + \int_{\partial\Omega_{2}} \frac{\partial P(x,y)}{\partial n}dS_{x} =-1.$ In the limit of two well separated windows ($|x_1-x_2| \gg 1)$, using the condition for the flux, we get for each window i=1,2, $\int_{\partial \Omega_{i}} \frac{\partial P(x,y))}{\partial n}dS_{x}= -\pi A_i$ (the minus sign is due to the outer normal orientation), thus $\pi A_1+\pi A_2=1.$ We finally obtain the system

$\frac{G(x_1,x_0)-G(x_2,x_0)}{\{\log|x_1-x_2|-\log \epsilon\}} +(A_2-A_1) =0$

$A_1+ A_2=\frac{1}{\pi}.$

Finally, the absorbing probabilities are given by $P_2=\pi A_2=\frac12+ \frac{\pi}{2}\frac{G(x_1,x_0)-G(x_2,x_0)}{\{\log|x_1-x_2|-\log \epsilon\}} =\frac12- \frac{1}{4}\frac{\ln \frac{|x_1-x_0|\left|x_1-\bar{x_0}\right|}{|x_2-x_0|\left|x_2-\bar{x_0}\right|}}{\{\log|x_1-x_2|-\log \epsilon\}}$and$P_1=\pi A_1=\frac12+ \frac{\pi}{2}\frac{G(x_2,x_0)-G(x_1,x_0)}{\{\log|x_1-x_2|-\log \epsilon\}} =\frac12- \frac{1}{4}\frac{\ln \frac{|x_2-x_0|\left|x_2-\bar{x_0}\right|}{|x_1-x_0|\left|x_1-\bar{x_0}\right|}}{\{\log|x_2-x_1|-\log \epsilon\}}.$

These probabilities precisely depend on the source position $x_0$ and the relative position of the two windows. When one of the splitting probabilities (either $P_1$ or $P_2$) is known and fixed in [0,1], recovering the position of the source requires inverting the previous equations. For $P_2=\alpha \in [0,1],$ the position $x_0$ lies on the curve

$S_{source}=\{ x_0 \hbox{ such that} \, \frac{|x_1-x_0|\left|x_1-\bar{x_0}\right|}{|x_2-x_0|\left|x_2-\bar{x_0}\right|}=\exp \left((4\alpha-2)\{\log|x_1-x_2|-\log \epsilon\} \right) \}.$ To conclude knowing the splitting probability between two windows is not enough to recover the exact the source $x_0$ , because it lies on one dimensional curve solution. However the direction can be obtained by simply checking which one of the two probability is the highest.

## Inverse probleme: computing source location by triangulation

To reconstruct the location of a source from the measured fluxes, at least three windows are needed. With two windows only, a source located on the line perpendicular to one of the connecting windows would, for example generate the same splitting probability $P_1=P_2$, leading to a one dimensional curve degeneracy for the reconstructed source positions $x_0$. Reconstructing the source location $x_0$ from the splitting probabilities,[7] requires to invert a system equation. Starting from the general solution as given above, we get for three windows:

$P(x,x_0)=G(x,x_0) +A_1 \{\log|x-x_1|-\log \epsilon\}+A_2 \{\log|x-x_2|-\log \epsilon\}+A_3 \{\log|x-x_3|-\log \epsilon\} +C,$

where $A_1,A_2,A_3,C$ are constants to be determined. The three absorbing boundary conditions for$P(x,x_0)$leads to the system of equations

$G(x_1,x_0) +A_2 \{\log|x_1-x_2|-\log \epsilon\}+A_3 \{\log|x_1-x_3|-\log \epsilon\} +C=0$

$G(x_2,x_0) +A_1 \{\log|x_2-x_1|-\log \epsilon\} +A_3 \{\log|x_2-x_3|-\log \epsilon\}+C=0$

$G(x_3,x_0) +A_1 \{\log|x_1-x_3|-\log \epsilon\}+A_2 \{\log|x_2-x_3|-\log \epsilon\} +C=0,$

with the normalization condition for the fluxes $\pi A_1+\pi A_2+\pi A_3=1.$ With $\Delta_{123}=\left(\log\frac{d_{13}d_{12}}{d_{32}\epsilon}\right)^2-4\log \frac{d_{12}}{\epsilon} \log \frac{d_{13}}{\epsilon}$, and using the notations $d_{ij}=|x_i-x_j|$, we obtain the solution

$A_1=\frac{1}{\pi}-A_2 -A_3.$

$A_2=\frac{ \log\frac{d_{13}d_{12}}{d_{32}\epsilon}(G_{30}-G_{10} +\frac{1}{\pi} \log\frac{d_{13}}{\epsilon})-( G_{1 0}-G_{2 0}+\frac{1}{\pi}\log\frac{d_{12}}{\epsilon})\log\frac{d_{13}^2}{\epsilon^2 })}{\Delta_{123}}.$

$A_3=\frac{ \log\frac{d_{13}d_{12}}{d_{32}\epsilon}(G_{20}-G_{10} +\frac{1}{\pi} \log\frac{d_{12}}{\epsilon})-( G_{1 0}-G_{3 0}+\frac{1}{\pi}\log\frac{d_{13}}{\epsilon}) \log\frac{d_{12}^2}{\epsilon^2})}{\Delta_{123}}$

To conclude, we can specify the flux values $\alpha>0$ and $\beta>0$ such that$\alpha+\beta<1$, with $\alpha=\pi A_1=\int_{\partial \Omega_{1}} \frac{\partial P(x,y))}{\partial n}dS_{x}$and $\beta=\pi A_2=\int_{\partial \Omega_{2}} \frac{\partial P(x,y))}{\partial n}dS_{x}.$ We now provide the final step which is to express the source $x_0$with respect to the flux coefficients $A_1,A_2,A_3.$With three windows in a x-y plane in the confirguration $x_1=(0, 0, 0), x_2=(d, 0, 0)\, \hbox{and} \, x_3=(e, f, 0)$ (window 1 is at the origin and window 2 is on the x-axis). Then, using the leading order from the expansion of the fluxes, the location of the source is the solution of the three non-linear equations

$\gamma_1^2 = (x_0^{(1)})^2 + (x_0^{(2)})^2 + (x_0^{(3)})^2$

$\gamma_2^2 = (d-x_0^{(1)})^2 + (x_0^{(2)})^2 + (x_0^{(3)})^2$

$\gamma_3^2 = (e-x_0^{(1)})^2 + (f-x_0^{(2)})^2 + (x_0^{(3)})^2,$

where $\gamma_i=\frac{2\epsilon}{\pi\Phi_i}$. Solving for the coordinates of \x_0 and requiring that x_0^{3}>0, leads to the analytical solution$x_0^{(1)}=\frac{d^2+\gamma_1-\gamma_2}{2d}$

$x_0^{(2)}= \frac{1}{2df}\left[d(e^2+f^2+\gamma_1-\gamma_3)-e(d^2+\gamma_1-\gamma_2)\right]$

$x_0^{(3)}= \frac{1}{2df}\left[(e^2+f^2)(\{\gamma_1-\gamma_2\}^2-d^4)+2de(e^2+f^2+\gamma_1-\gamma_3)(d^2+\gamma_1-\gamma_2)\right. \left.-d^2(e^4+f^4+[\gamma_1-\gamma_3]^2+2e^2[f^2+2\gamma_1-\gamma_2-\gamma_3]-2f^2[\gamma_2+\gamma_3])\right]^{1/2}.$


Thus the position can be recovered from the steady-state fluxes[8] [9] [10].

In general, with N windows, no explicit inverse is easily available, hence numerical procedures are used to find the position of the source x_0 at any order in $\epsilon$. To find the position of the source $x_0$ from the measured fluxes $\Phi_i, i=1...N,$we need to invert the previous Eqs. Each of these equations describes a non-planar surface $S_{i}$in three dimensions, corresponding to window i and intersecting the half-plane (in the case of the windows located on the half-plane) or the unit ball (in the case of the windows located on the ball). Each pair of surfaces $S_{i}$and $S_{j}$ intersect, forming three-dimensional curves $C_{ij}$and all of these curves intersect at the location of the source x_0. In the case of N>3 windows, the system is overdetermined and we can simply choose any combination k, l and m of three fluxes from the N available. Any combination leads to the same source position.

## References

1. Kolodkin, A. L.; Tessier-Lavigne, M. (2010-12-01). "Mechanisms and Molecules of Neuronal Wiring: A Primer". Cold Spring Harbor Perspectives in Biology. 3 (6): a001727–a001727. doi:10.1101/cshperspect.a001727. ISSN 1943-0264.
2. Blockus, Heike; Chédotal, Alain (August 2014). "The multifaceted roles of Slits and Robos in cortical circuits: from proliferation to axon guidance and neurological diseases". Current Opinion in Neurobiology. 27: 82–88. doi:10.1016/j.conb.2014.03.003. ISSN 0959-4388.
3. Dobramysl, U., & Holcman, D. (2018). Reconstructing the gradient source position from steady-state fluxes to small receptors. Scientific reports, 8(1), 1-8.
4. Shukron, O., Dobramysl, U., & Holcman, D. (2019). Chemical Reactions for Molecular and Cellular Biology. Chemical Kinetics: Beyond The Textbook, 353.
5. Pillay, S.; Ward, M. J.; Peirce, A.; Kolokolnikov, T. (January 2010). "An Asymptotic Analysis of the Mean First Passage Time for Narrow Escape Problems: Part I: Two-Dimensional Domains". Multiscale Modeling & Simulation. 8 (3): 803–835. doi:10.1137/090752511. ISSN 1540-3459.
6. John., Crank, (1956). The mathematics of diffusion. Clarendon. OCLC 1120831223.{{cite book}}: CS1 maint: extra punctuation (link) CS1 maint: multiple names: authors list (link)
7. 1937-, Schuss, Zeev, (2009). Theory and applications of stochastic processes : an analytical approach. Springer. ISBN 1-4614-2542-5. OCLC 1052813540.{{cite book}}: CS1 maint: extra punctuation (link) CS1 maint: multiple names: authors list (link) CS1 maint: numeric names: authors list (link)
8. Dobramysl, U., & Holcman, D. (2019). Triangulation sensing: how cells recover a source from diffusing particles in three dimensions. arXiv preprint arXiv:1911.02907.
9. Dobramysl, U., & Holcman, D. (2018). Mixed analytical-stochastic simulation method for the recovery of a Brownian gradient source from probability fluxes to small windows. Journal of computational physics, 355, 22-36.
10. DOBRAMYSL, Ulrich et HOLCMAN, David. Triangulation Sensing to Determine the Gradient Source from Diffusing Particles to Small Cell Receptors. Physical Review Letters, 2020, vol. 125, no 14, p. 148102.