Open Access
Add to BibSonomyAbstract
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
1. Introduction
The excited-state lifetime of a fluorescent 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 reflection by the surface [3]. Lifetime changes of emitters have been used for mapping the local density of optical states (LDOS) in the near field of cavities or nano-antennas [4,5]. In disordered scattering media, the fluorescence lifetime exhibits fluctuations due to changes of the local environment [6–11]. It has been demonstrated that even in a multiple scattering regime the amplitude of the fluctuations is sensitive to near-field 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 fluctuations are also influenced by the orientation of the transition dipole [13], leading to the development of new fluorescence polarization-based imaging techniques.
Near-field 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 amplification that arises from the inversion of evanescent waves since the noisy measurements are multiplied by exponentially increasing functions [14]. Second, in near-field imaging the measurements must be phase sensitive which is notoriously difficult 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 single-molecule sensitivity. Since the lifetime of a fluorophore 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 first 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-field 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.
2. Numerical model
In the weak coupling regime, the relative change in the decay rate Γ of a fluorescent source located at position r0 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, peg is the transition dipole between the excited and ground states, G⃡s is the scattered part of the dyadic Green function, and u = peg/|peg| is the normalized dipole direction. From this expression, it is apparent that decay rate data encode the flourophore position through the imaginary part of the scattered field Es(r0, ω) = μ0ω2G⃡s (r0, r0, ω)peg at its position.We consider here that the fluorophore 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 rj, 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]
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 field Es(r; r0, u, ω) at point r and frequency ω in terms of the exciting fields Eexc(rj; r0, u, ω) at each scatterer position. Both depend on the position r0, dipole direction u, and emision frequency ω of the fluorophore. The self-consistent system of M equations (3) gives the exciting fields at each scatterer position rj. The exciting field Eexc(rj; r0, u, ω) is equal to the sum of the incident field Einc(rj; r0, u, ω) at rj (radiated by the fluorophore with dipole direction u placed at r0) and the scattered fields at rj due to all scatterers at rm ≠ rj. 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 efficiently 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 fluorophore 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̂ is the tensor product of r̂ with itself, I⃡ is the unit tensor, and are the Hankel functions of the first 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 field 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 (kls)−1, where ℓs = 1/[ρσs(ω)] is the scattering mean free path. In this last expression, ρ denotes the density of scatterers and σs(ω) = (k3/8)|α(ω)|2 denotes the scattering cross section. In our experiments, ω0 = 2 · 106 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) = (1GHz, 2GHz), (kls)−1 varies from 0.23 to 0.18. These are frequencies close to the resonance frequency ω0. For frequencies in the range (δ1, δ2) = (3GHz, 5GHz), (kls)−1 decreases and varies from 0.11 to 0.09. For frequencies in the range (δ1, δ2) = (10GHz, 12GHz), (kls)−1 decreases even more and varies from 0.02 to 0.01.3. Numerical experiments
The lifetime of a fluorophore 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 configuration of scatterers but different fluorophore 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.
Fig. 1 Top row: Different fluorophore positions and orientations (r0, θ) 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 Γ(r0, θ, ω) of the fluorophores shown in the top row.
To compute the decay rate of a dipole source we have introduced the orientation parameter θ that defines the angle between the dipole direction u and the x-axis, that is, cos(θ) = u · x̂. For a fluorophore placed at r0 = (x, y) with orientation parameter θ, we have solved the M × M linear system (3) for the exciting fields, and we have computed the scattered field at r0 using Eq. (2). M=20 in all the numerical experiments shown in the paper. The fluorescent decay rates d(ω) = Γ(r0, θ, ω) 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 fluorophore and the environment is very sensitive to the frequency. This figure motivates the idea of recovering the fluorophore’s position r0 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 finding the position vector r0 and the dipole direction θ that satisfy the nonlinear equations
We consider a 2D configuration with r0 = (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 define a domain of interest (DOI) of size 1 × 1μm2 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 fluorophore using a standard microscope fluorescence 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 efficient gradient method that is easy to use. We note, though, that this is a local method that finds 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 finding 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.3.1. Noiseless data
In the first 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 first the best frequency range for which the solution to Eq. (6) is found more reliably. In Fig. 2, the scattering strength (kls)−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μm2 disordered medium in all the 6 experiments shown in the figure. The positions of the scatterers are shown by dots and the real fluorophore position by a star. The estimated position of the fluorophore 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 (kls)−1 = 0.18 − 0.23 (left column) the estimates are very precised, for (kls)−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.
Fig. 2 Reconstructions of the source positions in the same disordered medium as in Fig. 1. We consider a 1 × 1μm2 DOI (black square), and noiseless data from 3 frequencies. The scattering strength regimes are from left to right : (kls)−1 = 0.18 − 0.23, (kls)−1 = 0.09 − 0.11, and (kls)−1 = 0.01 − 0.02. The dipole source positions are shown by stars. The orientations are θ0 = 0 rad (top row) and θ0 = π/2 rad (bottom rows). The estimates of the source positions are shown by circles, and their orientations by arrows.
3.2. 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) = (1GHz, 2GHz), so the scattering regime is (kls)−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.
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 (kls)−1 = 0.18 − 0.23.
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.
Fig. 4 Reconstructions of the positions and orientations of dipole sources for different realizations of a disordered medium. We use data from 3 frequencies corrupted with 1% of additive noise. The scattering strength regime is (kls)−1 = 0.18 − 0.23.
3.3. 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 fluctuations around the true positions of the scatterers for the reconstructions.
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.
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 (kls)−1 = 0.18 − 0.23.
Figure 6 shows the results for the same examples as in Fig. 5 but using data from 9 frequencies within the range (δ1, δ2) = (1GHz, 2GHz), 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.
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.
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 when 3 frequencies (dot-dashed lines) and 9 frequencies (solid lines) are used. When 9 frequencies are used, we get a better defined local minimum of the residual. In the left plot of Fig. 7 we show the results for 3 frequencies (circle) and 9 frequencies (diamond). The true position of the fluorophore is depicted by a star.
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 (kls)−1 = 0.18 − 0.23.
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.
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.
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 fluorescent molecule in a strongly scattering medium. We have used the vector Foldy-Lax equations to calculate the interaction between a dipole fluorescent 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 find the values of five 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 fluorescence 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.
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. [CrossRef]
2. K. Drexhage, “Influence of a dielectric interface on fluorescence decay time,” J. Lumin. 1, 693–701 (1970). [CrossRef]
3. R.R. Chance, A. Prock, and R. Silbey, “Molecular fluorescence and energy transfer near interfaces,” Adv. Chem. Phys. 37, 1–65 (1978). [CrossRef]
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). [CrossRef] [PubMed]
5. M. Frimmer, Y. Chen, and A.F. Koenderink, “Scanning emitter lifetime imaging microscopy for spontaneous emission control,” Phys. Rev. Lett. 107, 123602 (2011). [CrossRef] [PubMed]
6. E.A. Donley and T. Plakhotnik, “Luminescence lifetimes of single molecules in disordered media,” J. Chem. Phys. 114, 9993–9997 (2001). [CrossRef]
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). [CrossRef] [PubMed]
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). [CrossRef]
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). [CrossRef] [PubMed]
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). [CrossRef]
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). [CrossRef] [PubMed]
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). [CrossRef]
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). [CrossRef]
14. G. Derveaux, G. Papanicolaou, and C. Tsogka, “Resolution and denoising in near-field imaging,” Inverse Problems 22, 1437–1456 (2006). [CrossRef]
15. A. Chai, M. Moscoso, and G. Papanicolaou, “Array imaging using intensity-only measurements,” Inverse Problems 27, 015005 (2011). [CrossRef]
16. N. Irishina, M. Moscoso, and R. Carminati, “Source location from fluorescence lifetime in disordered media,” Optics Letters , 37, 951–953 (2012). [CrossRef] [PubMed]
17. J.M. Wylie and J.E. Sipe, “Quantum electrodynamics near an interface,” Phys. Rev. A 30, 1185–1193 (1984). [CrossRef]
18. L.L. Foldy, “The multiple scattering of waves,” Phys. Rev. 67, 107–119 (1945). [CrossRef]
19. M. Lax, “Multiple scattering of waves,” Rev. Mod. Phys. 23, 287–310 (1951). [CrossRef]
20. S. Koc and W. C. Chew, “Calculation of acoustical scattering from a cluster of scatterers,” J. Acoust. Soc. Am. 103, 721–734 (1998). [CrossRef]
21. P. de Vries, D.V. van Coevordden, and A. Lagendijk, “Point scatterers for classical waves,” Rev. Mod. Phys. 70, 447–466 (1998). [CrossRef]
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).
