Recovering fluorophore location and orientation from lifetimes

In this paper, we study the possibility of using lifetime data to estimate the position and orientation of a fluorescent dipole source within a disordered medium. The vector Foldy-Lax equations are employed to calculate the interaction between the fluorescent source and the scatterers that are modeled as point-scatterers. The numerical experiments demonstrate that if good prior knowledge about the positions of the scatterers is available, the position and orientation of the dipole source can be retrieved from its lifetime data with precision. If there is uncertainty about the positions of the scatterers, the dipole source position can be estimated within the same level of uncertainty. © 2013 Optical Society of America OCIS codes: (290.3200) Inverse scattering; (260.2510) Fluorescence. References and links 1. See K. Suhling, P.W. French and D. Philipps, “Time-resolved fluorescence microscopy,” Photochem. Photobiol. Sci. 4, 13–22 (2005) and references therein. 2. K. Drexhage, “Influence of a dielectric interface on fluorescence decay time,” J. Lumin. 1, 693-701 (1970). 3. R.R. Chance, A. Prock and R. Silbey, “Molecular fluorescence and energy transfer near interfaces,” Adv. Chem. Phys. 37, 1–65 (1978). 4. J.P. Hoogenboom, G. Sanchez-Mosteiro, G. Colas des Francs, D. Heinis, G. Legay, A. Dereux and N.F. van Hulst, “The single molecule probe: nanoscale vectorial mapping of photonic mode density in a metal nanocavity,” Nano Lett. 9, 1189–1195 (2009). 5. M. Frimmer, Y. Chen and A.F. Koenderink, “Scanning emitter lifetime imaging microscopy for spontaneous emission control,” Phys. Rev. Lett. 107, 123602 (2011). 6. E.A. Donley and T. Plakhotnik, “Luminescence lifetimes of single molecules in disordered media,” J. Chem. Phys. 114, 9993–9997 (2001). 7. R.A.L. Vallée, N. Tomczak, L. Kuipers, G.J. Vancso and N.F. van Hulst, “Single molecule lifetime fluctuations reveal segmental dynamics in polymers,”Phys. Rev. Lett. 91, 038301 (2003). 8. L.S. Froufe-Pérez, R. Carminati and J.J. Sáenz, “Fluorescence decay rate statistics of a single molecule in a disordered cluster of nanoparticles,” Phys. Rev. A 76, 013835 (2007). 9. M. D. Birowosuto, S. E. Skipetrov, W. L. Vos and A. P. Mosk, “Observation of spatial fluctuations of the local density of states in random media,” Phys. Rev. Lett. 105, 013904 (2010). 10. V. Krachmalnicoff, E. Castanié, Y. De Wilde and R. Carminati, “Fluctuations of the local density of states probe localized surface plasmons on disordered metal films,” Phys. Rev. Lett. 105, 183901 (2010). 11. R. Sapienza, P. Bondareff. R. Pierrat, B. Habert, R. Carminati and N. F. van Hulst, “Long-Tail statistics of the Purcell factor in disordered media driven by near-field interactions,” Phys. Rev. Lett. 106, 163902 (2011). 12. A. Cazé, R. Pierrat and R. Carminati, “Near-field interactions and nonuniversality in speckle patterns produced by a point source in a disordered medium,” Phys. Rev. A 82, 043823 (2010). 13. L.S. Froufe-Pérez and R. Carminati, “Lifetime fluctuations of a single emitter in a disordered nanoscopic system: The influence of the transition dipole orientation,” Phys. Stat. Sol. (a) 205, 1258–1265 (2008). 14. G. Derveaux, G. Papanicolaou and C. Tsogka, “Resolution and denoising in near-field imaging,” Inverse Problems 22, 1437-1456 (2006). #176443 $15.00 USD Received 18 Sep 2012; revised 15 Nov 2012; accepted 16 Nov 2012; published 4 Jan 2013 (C) 2013 OSA 14 January 2013 / Vol. 21, No. 1 / OPTICS EXPRESS 421 15. A. Chai, M. Moscoso, and G. Papanicolaou, “Array imaging using intensity-only measurements,” Inverse Problems 27 , 015005 (2011). 16. N. Irishina, M.Moscoso and R. Carminati, “Source location from fluorescence lifetime in disordered media,” Optics Letters, 37, 951–953 (2012). 17. J.M. Wylie and J.E. Sipe, “Quantum electrodynamics near an interface,” Phys. Rev. A 30, 1185–1193 (1984). 18. L.L. Foldy, “The multiple scattering of waves,” Phys. Rev. 67, 107–119 (1945). 19. M. Lax, “Multiple scattering of waves,” Rev. Mod. Phys. 23, 287–310 (1951). 20. S. Koc and W. C. Chew, “Calculation of acoustical scattering from a cluster of scatterers,” J. Acoust. Soc. Am. 103, 721–734 (1998). 21. P. de Vries, D.V. van Coevordden and A. Lagendijk, “Point scatterers for classical waves,” Rev. Mod. Phys. 70, 447–466 (1998). 22. A. Ruszczynski, Nonlinear Optimization (Princeton University Press, Princeton, 2006). 23. E. R. Hansen, Global Optimization using Interval Analysis (Marcel Dekker, New York, 1992). 24. R. Horst, P. M. Pardalos, and N. V. Thoai, eds., Introduction to Global Optimization (Kluwer Academic, 2nd ed., 2000).


Introduction
The excited-state lifetime of a fluorescen emitter is a sensitive probe of its chemical and electromagnetic environment.In biology, this property has led to the development of sensitive imaging techniques [1].In nanophotonics, the change in the lifetime of emitters close to metallic surfaces has been demonstrated long ago [2], and explained in terms of classical electromagnetic interaction between a dipole source and its reflectio by the surface [3].Lifetime changes of emitters have been used for mapping the local density of optical states (LDOS) in the near fiel of cavities or nano-antennas [4,5].In disordered scattering media, the fluorescenc lifetime exhibits fluctuation due to changes of the local environment [6][7][8][9][10][11].It has been demonstrated that even in a multiple scattering regime the amplitude of the fluctuation is sensitive to near-fiel interactions between the emitter and the nearby scatterers, thus providing a way to get structural information at subwavelength scale in complex media [11,12].The lifetime fluctuation are also influence by the orientation of the transition dipole [13], leading to the development of new fluorescenc polarization-based imaging techniques.
Near-fiel imaging techniques for the detection of objects are receiving great attention.It should be emphasized, however, that these techniques suffer from two drawbacks.First, to achieve good resolutions they should include high spatial frequency evanescent components of the waves, and these exponentially small components have to be measured very close to the object.Furthermore, if evanescent components are used, the problem is noise amplificatio that arises from the inversion of evanescent waves since the noisy measurements are multiplied by exponentially increasing functions [14].Second, in near-fiel imaging the measurements must be phase sensitive which is notoriously difficul in optics.To overcome this problem, new imaging methods from measurements only of the power received at the detectors have been developed [15].The authors in [15] observed that the power is a linear function of a rank one matrix associated with the unknown source positions and, thus, the reconstruction is formulated as a rank minimization problem.
Alternatively, in this paper we propose a novel strategy based on lifetime data.Lifetime measurements provide a novel, intensity-independent mechanism in imaging with singlemolecule sensitivity.Since the lifetime of a fluorophor in a scattering environment strongly depends on its position and orientation, it should be possible to recover them from the lifetime data, provided that a certain degree of prior information is available (in particular the geometry of the environment).The purpose of this work is to study this question on the basis of numerical experiments.
A firs proof of principle has been given recently in the case of the scalar approximation to the problem, showing the feasibility of the source location in a medium consisting of discrete point-like scatterers [16].Here, we address the problem in a much more general context.First, we perform numerical experiments taking into account the polarization of the light field This is a key point since near-fiel interactions strongly depend on the orientation of the dipole source (emitter) with respect to nearby scatterers.Second, we study the robustness of this method by reducing the degree of knowledge on the geometry of the scattering medium.To this end, we introduce some uncertainty on the positions of the scatterers and we quantify the impact on the accuracy of the results.
The paper is organized as follows.In section 2, we introduce the numerical model.In section 3, we present numerical experiments with the reconstructions of the positions and orientations of a dipole source in different environments.We consider different scattering strength regimes and we relax the prior knowledge on the positions of scatterers.Section 4 contains our conclusions.

Numerical model
In the weak coupling regime, the relative change in the decay rate Γ of a fluorescen source located at position r 0 due to the interaction with the surrounding medium is given by [17] In Eq. ( 1), Γ 0 is the decay rate in free space, ω is the emission frequency, p eg is the transition dipole between the excited and ground states, ↔ G s is the scattered part of the dyadic Green function, and u = p eg /|p eg | is the normalized dipole direction.From this expression, it is apparent that decay rate data encode the flourophor position through the imaginary part of the scattered fiel E s (r 0 , ω) = µ 0 ω 2 ↔ G s (r 0 , r 0 , ω)p eg at its position.We consider here that the fluorophor is immersed within a discrete multiple scattering medium that is composed of M scatterers with polarizability α(ω) and that are small compared to the wavelength.The scatterers are distributed randomly in the free space an have positions r j , j = 1, . . ., M. Wave propagation and multiple scattering in such a medium can be formulated in terms of the so-called Foldy-Lax equations [18,19] E exc (r j ; r 0 , u, ω) = E inc (r j ; r 0 , u, ω) Here, ↔ G 0 is the dyadic free space Green function, k = ω/c is the wavenumber, and c is the speed of light in free space.The summations are extended over the discrete set of scatterers.Eq. ( 2) gives the scattered fiel E s (r; r 0 , u, ω) at point r and frequency ω in terms of the exciting field E exc (r j ; r 0 , u, ω) at each scatterer position.Both depend on the position r 0 , dipole direction u, and emision frequency ω of the fluorophore The self-consistent system of M equations (3) gives the exciting field at each scatterer position r j .The exciting fiel E exc (r j ; r 0 , u, ω) is equal to the sum of the incident fiel E inc (r j ; r 0 , u, ω) at r j (radiated by the fluorophor with dipole direction u placed at r 0 ) and the scattered field at r j due to all scatterers at r m = r j .This system can be written in matrix form and solved numerically by standard techniques if the number of scatterers M is not too large.In our experiments, M is small (M = 20) and we solve the linear system of equations (3) with a standard LU factorization technique.For M large, Eqs.(3) can be solved efficientl using, for example, the Fast Multipole Method [20].
In the numerical experiments shown in the next section we will consider a 2D configuration The scatterers lie in the xy plane, and the dipole moment of the fluorophor is parallel to this plane (p polarization).In this case, the free-space dyadic Green function is given by Here, R = ||r − r'||, r = r−r' ||r−r'|| , r ⊗ r is the tensor product of r with itself, ↔ I is the unit tensor, and H (1) i are the Hankel functions of the firs type and order i = 0, 1.The scatterers are described by the resonant polarizability where ω 0 is the resonance frequency and γ is the linewidth.This form of the polarizability describes a lossless scatterer, and includes self interaction.This is the reason why only the exciting fiel on each scatterer is involved in Eq. ( 3).Including self-interaction (radiation reaction) is also a requirement for energy conservation, or equivalently the optical theorem [21].We note that multiple scattering regimes with different strengths can be explored by adjusting the detuning frequency δ = ω − ω 0 , with δ ω 0 .As a measure of the scattering strength, we use (kl s ) −1 , where s = 1/[ρσ s (ω)] is the scattering mean free path.In this last expression, ρ denotes the density of scatterers and σ s (ω) = (k 3 /8)|α(ω)| 2 denotes the scattering cross section.In our experiments, ω 0 = 2 • 10 6 GHz (λ 0 = 0.94 µm) and γ = 5 GHz.In order to study scattering regimes with different strengths, we select appropriate ranges of frequencies.For frequencies ω in the range (δ 1 , δ 2 ) = (1 GHz, 2 GHz), (kl s ) −1 varies from 0.23 to 0.18.These are frequencies close to the resonance frequency ω 0 .For frequencies in the range (δ 1 , δ 2 ) = (3 GHz, 5 GHz), (kl s ) −1 decreases and varies from 0.11 to 0.09.For frequencies in the range (δ 1 , δ 2 ) = (10 GHz, 12 GHz), (kl s ) −1 decreases even more and varies from 0.02 to 0.01.

Numerical experiments
The lifetime of a fluorophor in a disordered medium strongly depends on its position and orientation with respect to the scatterers [6,8,13].To illustrate this phenomenon, we have carried out a set of numerical experiments for the same configuratio of scatterers but different fluorophor positions and orientations.The top row in Fig. 1 shows with dots the positions of the scatterers, and with a star the position of the fluorophore The fluorophore s orientations are indicated with arrows.
To compute the decay rate of a dipole source we have introduced the orientation parameter θ that define the angle between the dipole direction u and the x-axis, that is, cos(θ ) = u • x.For a fluorophor placed at r 0 = (x, y) with orientation parameter θ , we have solved the M × M linear system (3) for the exciting fields and we have computed the scattered fiel at r 0 using Eq.(2).M=20 in all the numerical experiments shown in the paper.The fluorescen decay rates d(ω) = Γ(r 0 , θ , ω) at different frequencies ω are then computed using Eq. ( 1).The results are shown in the bottom row of Fig. 1.The decay rate strongly depends on the frequency ω because the interaction between the fluorophor and the environment is very sensitive to the frequency.This figur motivates the idea of recovering the fluorophore s position r 0 and its orientation θ from the knowledge of the decay rate at different frequencies.In the next numerical experiments we show that this can, in fact, be accomplished easily.Mathematically, we formulate the problem as that of findin the position vector r 0 and the dipole direction θ that satisfy the nonlinear equations We consider a 2D configuratio with r 0 = (x, y) and, thus, we only seek for the values of three unknowns: (x, y) and θ .Hence, we only need, in principle, data corresponding to three different frequencies to recover these three unknowns by solving Eq. ( 6).If desired, we can use data from more frequencies and compute the nonlinear least-squares solution to Eq. ( 6) for estimates of those unknowns.This will be the case when the data is contaminated with noise and/or there is uncertainty in the position of the scatterers.
In the numerical experiments shown below, we defin a domain of interest (DOI) of size 1×1 µm 2 around the true source position to restrict the search area.This can be done in practice using, for example, a coarse-grain for locating the fluorophor using a standard microscope fluorescenc intensity image.Once the DOI is defined the position and orientation of the dipole source are found by minimizing the residual using a standard iterative gradient method.We use the Matlab Optimization Toolbox that contains an efficien gradient method that is easy to use.We note, though, that this is a local method that find only the nearest local minimum to the initial guess.Other local methods for solving nonlinear equations can be used as well [22].A good alternative for findin the global minimum, at the cost of a much higher computational cost, is the use of global methods such as, for example, interval methods, simulating annealing or Monte-Carlo methods [23,24].These methods can only be used for small problems, though.

Noiseless data
In the firs example we assume perfect data and exact knowledge of the positions of the M = 20 scatterers, as well as their polarizability.We do not know neither the position nor the dipole orientation of the fluorophore Since different ranges of frequency correspond to distinct regimes of interaction between the source and its environment, we investigate firs the best frequency range for which the solution to Eq. ( 6) is found more reliably.In Fig. 2, the scattering strength (kl s ) −1 decreases from the left column to the right column.For each frequency range (or each scattering regime), we use synthetic data corresponding to only N ω = 3 frequencies within that range of frequency.We keep the same 3.7 × 3.7 µm 2 disordered medium in all the 6 experiments shown in the figure The positions of the scatterers are shown by dots and the real fluorophor position by a star.The estimated position of the fluorophor found by solving Eq. ( 6) is shown by a circle.The orientations of both, the real and reconstructed dipole source, are represented by arrows.As reported in [16], we observe that the realiability of the reconstruction deteriorates as the scattering strength decreases.While for (kl s ) −1 = 0.18 − 0.23 (left column) the estimates are very precised, for (kl s ) −1 = 0.01−0.02 the estimates are not as good.This is so because lifetime data is more sensitive to small changes in the position or orientation of the dipole source when its interaction with the environment is stronger and, hence, the global minimum of Eq. ( 7) is better defined When the interaction between the dipole source an the environment is weak, many positions and orientations of the dipole source give rise to nearly the same lifetime data.

Noisy data
In the next examples we investigate the robusteness of the proposed method to noise in the data.We still assume perfect knowledge of the positions of the scatterers and we use 3 frequencies within the range (δ 1 , δ 2 ) = (1 GHz, 2 GHz), so the scattering regime is (kl s ) −1 = 0.18 − 0.23.
Figure 3 shows the results when the data is corrupted with 1% (left column), 3% (middle column) and 8% (right column) of additive gaussian noise.The positions and orientations of the sources are the same as in Fig. 2. The positions of the scatterers are represented by dots, and the positions of the real and recovered sources by a star and a circle, respectively.Their orientations are represented by arrows.As expected, the estimates of the source positions and orientations deteriorates with noise.We observe quite good estimates for noise below 3%, though.
Figure 4 displays results corresponding to more numerical experiments.We now consider different realizations of the disordered media.In all the cases, 1% of noise is added to the data.We observe that the estimates of the positions and orientations of the dipole sources are always quite good.

Uncertainty with respect to the scatterers positions
In the previous examples we have assumed that the positions of the scatterers were known exactly.Under this assumption, the dipole source is localized and characterized correctly from lifetime data if the noise in the data is not too large.Next, we investigate the robustness of this method when the positions of the scatterers are not known exactly.We model this uncertainty by adding arbitrary random fluctuation around the true positions of the scatterers for the reconstructions.
x (µ m) y (µ m) Fig. 6.Same as Fig. 5, but using data from 9 frequencies.There is a ±λ /8 uncertainty on the scatterers positions and the data is corrupted with 1% of additive noise.In Fig. 5 we have assumed an uncertainty of ±λ /8 in the scatterer's positions.We still use data from 3 frequencies corrupted with 1% of gaussian noise.We now observe that in almost all the cases the estimates of the positions and orientations of the dipole sources are not satisfactory.
Figure 6 shows the results for the same examples as in Fig. 5 but using data from 9 frequencies within the range (δ 1 , δ 2 ) = (1 GHz, 2 GHz), instead of only 3 frequencies.When data from more frequencies are used, the positions and orientations are better estimated.We note, however, that beyond a certain number of frequencies the values of the estimates do not improve anymore.In general, when enough frequencies are used the distances between the true and estimated positions are of the same order as the uncertainty on the scatterer's positions.
To better understand the role of the number of frequencies on the quality of the results, we show in the middle and right plots of Fig. 7 the x and y cross-sections of the residual  Finally, in Fig. 8 we display the reconstructions of the same dipole sources shown in Fig. 5 but with an uncertainty of ±λ /4 in the scatterers positions.We use lifetime data from 9 frequencies corrupted with 1% of gaussian noise.We have applied now a global optimization method since multiply local minima in the cost functional exist within the DOI and, therefore, a gradient method fails.The orientation of the dipoles are well estimated in all the cases, and their positions are recovered within a precision of ±λ /4.

Conclusion
We have presented a novel intensity-independent strategy based on lifetime data that is able to localize and characterize with precision the emission of a fluorescen molecule in a strongly scattering medium.We have used the vector Foldy-Lax equations to calculate the interaction between a dipole fluorescen source and the surrounding scatterers in a 2D configuration The more general 3D case can be carried out, in principle, without modification The major difference is that in 3D we would have to fin the values of f ve unknowns (three spatial variables (x, y, z) and two directional variables (θ , φ )) and, therefore, the resolution of the nonlinear optimization problem is more complicated.The numerical experiments show that the position and orientation of a dipole source inside a predetermined domain of interest can be retrieved with high accuracy from multi-frequency fluorescenc lifetime data (except for the ambiguity under dipole inversion).The proposed strategy requires prior knowledge about the scatterers position.Our results indicate, however, that the precision of the recovered source position deteriorates only at the same rate as the uncertainty on the scatterer's position increases.

Fig. 1 .
Fig. 1.Top row: Different fluorophor positions and orientations (r 0 , θ ) within the same realization of a disordered medium.The fluorophore s orientations are indicated with an arrow.Bottom row: Frequency dependence of the decay rates Γ(r 0 , θ , ω) of the fluorophore shown in the top row.

Fig. 3 .Fig. 4 .
Fig. 3. Same as Fig. 2 but with data corrupted by noise.The noise levels (from left to right) are: 1%, 3% and 8%.The estimates of the source positions are shown by circles, and their orientations by arrows.The scattering strength regime is (kl s ) −1 = 0.18 − 0.23.

Fig. 5 .
Fig.5.Reconstructions of the same dipole sources as in Fig.4, but with an uncertainty on the scatterers positions of ±λ /8.We use data from 3 frequencies corrupted by 1% of noise.The scattering strength is (kl s ) −1 = 0.18 − 0.23.

Fig. 7 .
Fig. 7. Left panel: Reconstructions of a dipole source (shown by star) using 3 (circle) and 9 (diamond) frequencies.Middle and right panels: Residuals along the x and y axes, respectively.Dashed-dotted lines correspond to 3 frequencies, solid lines correspond to 9 frequencies.The data is corrupted with 1% of gaussian noise, and there is an uncertainty of ±λ /8 on the scatterers positions.The scattering strength is (kl s ) −1 = 0.18 − 0.23.
2 when 3 frequencies (dot-dashed lines) and 9 frequencies (solid lines) are used.When 9 frequencies are used, we get a better define local minimum of the residual.In the left plot of Fig.7we show the results for 3 frequencies (circle) and 9 frequencies (diamond).The true position of the fluorophor is depicted by a star.

Fig. 8 .
Fig. 8. Reconstructions of the same dipole sources as in Figs. 4 and 5, but using data from 9 frequencies (with 1% noise) and an uncertain on the scatterers positions of ±λ /4.