Apparent self-heating of individual upconverting nanoparticle thermometers

Individual luminescent nanoparticles enable thermometry with sub-diffraction limited spatial resolution, but potential self-heating effects from high single-particle excitation intensities remain largely uninvestigated because thermal models predict negligible self-heating. Here, we report that the common “ratiometric” thermometry signal of individual NaYF4:Yb3+,Er3+ nanoparticles unexpectedly increases with excitation intensity, implying a temperature rise over 50 K if interpreted as thermal. Luminescence lifetime thermometry, which we demonstrate for the first time using individual NaYF4:Yb3+,Er3+ nanoparticles, indicates a similar temperature rise. To resolve this apparent contradiction between model and experiment, we systematically vary the nanoparticle’s thermal environment: the substrate thermal conductivity, nanoparticle-substrate contact resistance, and nanoparticle size. The apparent self-heating remains unchanged, demonstrating that this effect is an artifact, not a real temperature rise. Using rate equation modeling, we show that this artifact results from increased radiative and non-radiative relaxation from higher-lying Er3+ energy levels. This study has important implications for single-particle thermometry.


Supplementary Note 1: Estimating R internal , R substrate , and R contact
In our thermal model, we assume that the entire nanoparticle is at a uniform temperature T NP . To validate this "lumped" assumption, we estimate an effective internal thermal resistance for the nanoparticle and compare this value to the largest external resistor, R air .
While the thermal resistance concept technically is not applicable in this case because the internal energy generation is volumetrically distributed throughout the particle, rather than occurring at a single point location, by approximating the nanoparticle as a sphere of diameter D with spatially uniform heat generation, we can define a suitable effective resistor with the correct K W -1 units -namely, the difference between the central and surface temperatures divided by the heat loading, which gives !"#$%"&' = 1 (4 !"#$%&'( ). We take a conservatively low value of k particle ≈ 1 W/m-K for the nanoparticle thermal conductivity 1 . Thus, for a 50 x 50 x 50 nm 3 particle, R int ≈ 2 x 10 6 K W -1 . Therefore, the Biot number Bi = R internal /R air is ~10 -3 << 1 and internal temperature gradients can be safely neglected to an excellent approximation.
R substrate is estimated as 1 2 !"#!$%&$' by approximating the particle as a heated disk on a semi-infinite medium, where k substrate is the substrate thermal conductivity. Because the majority of our experiments employ faceted 50 x 50 x 50 nm 3 nanorods with their hexagonal faces laying flat on a substrate, we take D to be 50 nm. Here, we use the bulk value of k glass ≈ 1 W m -1 K -1 and ignore possible phonon size effects, because we expect the average diameter of the contact area to be larger than the ~1 nm average phonon mean free path for glass at room temperature 2 . Under these approximations, we calculate R substrate ≈ 1 x 10 7 K W -1 for a 50 x 50 x 50 nm 3 particle on a borosilicate glass substrate (the primary experimental configuration in this work). Thus although we neglect R substrate for the purpose of our conservative thermal estimate, in reality the substrate may act as a significant heat sink, which if accounted for in this analysis can only make the particle colder.
R contact is the most challenging resistor to estimate due to the unknown contact area and limited theoretical understanding of thermal transport at nanoparticle-substrate interfaces.
Calculated values for a 1 µm diameter Si nanowire on Si substrate surrounded by air 3 range from ~10 4 -10 8 K W -1 . Similarly, room temperature experimental values for carbon nanofibers and multiwalled carbon nanotubes on various substrates 4,5 are between ~10 4 -10 7 K W -1 . If R contact in our experiments is on the high end of the range of reported values for these nanostructures with similar characteristic lengths, then R contact is comparable to R air and a significant portion of the heat will be dissipated through the air. A very large R contact would not change the estimated thermal resistance or temperature rise, since in this case R substrate + R contact >> R air holds true. The magnitude of R contact does impact the anticipated effects of changing other components of the thermal circuit. For example, if R contact is the limiting resistance, then we expect the apparent temperature rise to be insensitive to changes in R substrate , as discussed in the main text. If R contact is instead on the lower end of the range, a significant portion of the heat may flow through the contact. If R substrate is also small compared to R air (i.e. now R air >> R substrate + R contact ), R model and thus the estimated temperature rise can only decrease.

Supplementary Note 2: Single particle identification and characterization
To identify individual nanoparticles with characteristic dimensions far smaller than the diffraction-limited laser beam diameter, we followed a well-established statistical approach 6,7 in order to determine the characteristic brightness of a single particle. We obtained an APD scan that contains a total of approximately 300 luminescent emission spots. Supplementary Figure 1 shows a histogram of this data for the 50 x 50 x 50 nm 3 particles and the inset shows a representative portion of the APD scan. The characteristic brightness of a single particle can clearly be identified as ~6.5 kcnts per s. Given this information, when a new sample from the same batch of particles is imaged under the same conditions, we can rapidly identify single particles to use for subsequent thermal and intensity-dependent measurements. Supplementary Figure 1 also shows that the vast majority of the luminescent spots correspond to single nanoparticles, with a very small number of higher-intensity spots (> 10 kcnts/s) representing nanoparticle clusters, indicating that nanoparticle aggregation during the sample preparation process is minimal.  The green spectral emission of NaYF 4 :Yb 3+ ,Er 3+ UCNPs spans a wavelength range of approximately 515-565 nm. For the purpose of ratiometric thermometry, this emission band is typically separated into two sub-bands, a high-energy band of 515 nm < λ < 535 nm and a low energy band of 535 nm < λ < 560 nm. The high-energy band represents emission due to the 2 H 11/2 to 4 I 15/2 transition and the low-energy band represents emission due to the 4 F 9/2 to 4 I 15/2 transition. As noted in the main text, we exclude the peak observed at approximately 556 nm because of its known non-thermal origin. If instead that 556 peak is included in the second, low-energy band, the luminescence intensity ratio r displays a non-monotonic dependence on excitation intensity. As shown in

Supplementary Note 4: Error bars for ratiometric apparent temperature rise data
The error bars for all ratiometric apparent temperature rise data (i.e. the five different slopes from fitting ln(r) vs. 1/T (a linear relationship obtained by taking natural log of Eq. (1) in the main text) for each particle shown in Fig. 2(a). Because we are interested in calculating apparent temperature rises (as opposed to absolute temperatures), it is the variation in these slopes that is important for our uncertainty analysis. Therefore, for every measured r value, the change in the ratio with respect to the zero-power value, i.e. ∆r = r -r intrinsic (T = 296 K), where r intrinsic (T = 296 K) ≈ 0.289 from Fig. 2(c), is calculated. Two consecutive spectra were taken at each I exc , resulting in two ∆r values for every I exc . For each ∆r, five different values of the apparent temperature rise are calculated using each of the five fitted slopes. This approach results in in ten apparent temperature rise values calculated at each I exc . Finally, we plot the mean of these ten calculated values, and the error bars represent the standard deviation. Thus the error bars reflect both the variation between consecutive measurements and among particles in the same batch, which are the dominant sources of uncertainty in our measurements.

Supplementary Note 5: Error bars for luminescence lifetime apparent temperature rise data
The error bars for the luminescence lifetime apparent temperature rise data ( Fig. 3(b) and (d)) are calculated in the same manner as the error bars for the ratiometric apparent temperature rise data, except that the underlying model is now a linear fit of τ lum vs. T.
For every measured τ lum value, the change in the lifetime with respect to the zero-power value is calculated as ∆τ lum = τ lumτ lum,intrinsic (T = 296 K), where τ lum,intrinsic (T = 296 K) is obtained by extrapolating the τ lum (I exc ,T) data at T = 296 K. We consider the variation in the slopes of the lines fitted to the lifetime vs. temperature data for the five particles shown in Fig. 3(a) and (c). Here, three consecutive measurements were performed at every I exc , resulting in fifteen apparent temperature rise values for each I exc . Again, the mean of the fifteen calculated values is plotted, and the error bars represent the standard deviation.

Supplementary Note 6: Rate Equation Modeling
Differential rate equations (DREs) were used to model steady-state changes in population n i for each excited state i in the Yb 3+ and Er 3+ dopants. Our computational model 11,12 solves systems of coupled DREs that describe the rate at which each lanthanide 4f N manifold i is populated and depopulated by photon absorption, luminescence, energy transfer, and multiphonon relaxation. Supplementary Figure 4(a) displays the shift in population to higher-lying Er 3+ energy levels as the excitation intensity is increased, and Supplementary Fig. 4(b) shows the resulting spectral emission changes with excitation intensity. As indicated by the magenta dashed lines in Supplementary Fig. 4(b), the wavelength bands used to calculate the modeled luminescence intensity ratio are 515-535 nm and 535-555 nm. The wavelength cutoff for the second band is slightly higher than the corresponding experimental value of 548 nm to account for the fact that the experimental spectral peaks have finite widths, in contrast with the simulations. Thus, although the experimental emission spectra are integrated only to 548 nm, some of this integrated emission intensity comes from transitions centered at slightly longer wavelengths. Supplementary Fig. 4(c) quantifies the increase in emission originating from the population change in the 2 H 11/2 manifold relative to the 4 S 3/2 manifold as I exc is increased. Supplementary Fig. 4(d)  Here we present the derivation of the modeled curves plotted in Fig. 6(c)

of the main text
The steady-periodic solution to the governing differential equation for the frequencydependent temperature rise (main text Eq. (4)) can be obtained analytically using several different approaches. Here, we present a piecewise solution that takes advantage of the periodicity of Q(t). We start by defining the period of excitation, P, as P = 1/f exc . For times t such that NP < t < NP + P/2, where N = 0, 1, 2..., Q(t) has a constant value of Q max . For times t such that NP + P/2 < t < (N+1)P, Q(t) has a constant value of 0.
Implicitly, this analysis assumes that the rise and fall time of Q(t) is much faster than τ thermal or τ lum . The governing equation can therefore be written as: Both parts of Supplementary Equation (1) can be solved analytically using standard approaches. Since θ(t) must be continuous and has period P, the solutions to both parts of Supplementary Equation (1) must equal each other at t = NP + P/2, and also at t = (N+1)P.
Because we are interested in the steady-periodic solution, N does not matter and is henceforth set to zero. The resulting piecewise solution to Supplementary Equation (1) can be expressed as follows: An analogous governing equation can be written for the luminescence response. I lum increases with I exc , which is here a square wave with period P. The response of I lum to a step change in I exc is taken to be a first-order exponential relaxation, with a time constant τ lum . The governing equation for I lum (t) can thus be written as follows: By applying analogous stitching conditions every half period, the solution to Supplementary Equation (3) can be expressed as: In this analytical framework, we assume that both τ lum and I lum,max are constant, although in reality we observe that these quantities are modest functions of temperature. Figure   2(a) shows that the τ lum decreases from approximately 260 µs to 180 µs between 296 K and 400 K, or roughly 3000 ppm K -1 . Accounting for this effect induces a negligible shift in the characteristic frequency at which the apparent temperature rise begins to drop off as compared to the curves shown in Fig. 6(c). We also observe from our experiments that I lum,max decreases as a function of temperature. At steady state, this is equivalent to approximating !"# = • !"# , where we take γ to be independent of T. In reality, γ depends weakly on T, but this effect (also ~3000 ppm K -1 ) is negligible compared to the much larger contrast between I exc in the on and off states and thus has a similarly minimal impact on the modeled curves shown in Fig. 6(c). Consequently, we neglect these two second-order effects so that we can obtain an analytical solution for θ measured (f exc ).
In the main text, we qualitatively describe the physical phenomena that lead to a loss of sensitivity to R thermal when τ lum >> τ thermal . Here, we provide a more detailed, quantitative explanation of their effect on θ measured . The measurement is easiest to understand when τ lum à 0, such that luminescence is only emitted during the half cycle when the excitation laser is on, at an approximately constant value I lum,on . Again, I lum (t) is approximately a square wave between 0 and I lum,max , synchronized with the I exc (t) square wave between 0 and I exc,on. Eq. (5) in the main text then simplifies to As discussed in the main text, here !"#$%&"',! !"# →! is an average of the true θ(t) but averaged only over the half cycle when I exc is on. If τ thermal << 1/f exc , this will simply yield the DC value, θ DC . This θ DC is an important reference value for other regimes. For example, if τ thermal >> 1/f exc >> τ lum , the temperature cannot keep up with the laser heating oscillations and instead θ(t) à 0.5θ DC , so θ measured = 0.5θ DC (Fig. 6(b)). Thus, for τ lum à 0, the frequency at which θ measured transitions from θ DC to 0.5θ DC corresponds to f exc τ thermal ~1.
However, in the actual experimental regimes of greatest interest, we find that τ lum cannot simply be set to 0. Indeed, the most significant consequence of Eq. (5) occurs when τ lum >> 1/f exc >> τ thermal . In this regime, the spectrometer spectral ratio will also correspond to θ measured ≈ 0.5θ DC , even though in reality θ(t) remains a square wave between 0 and DC θ .
Due to the large τ lum , I lum (t) now remains essentially constant throughout the cycle, and Eq. (5) simplifies to In contrast to Supplementary Equation (5), now !"#$%&"',! !"# →! represents an average of the true θ(t) over the full cycle. Thus, even though the true θ in Fig. 6(a) still reaches θ DC , the spectrometer ratio will correspond to 0.5θ DC , but for a reason that has nothing to do with f exc τ thermal ~1. Here, this transition instead occurs when f exc τ lum ~1. As noted in the main text, the key takeaway is that there are two very different mechanisms that can cause a transition from θ measured = θ DC à 0.5θ DC with increasing f exc , resulting in a loss of sensitivity to R thermal when τ lum >> τ thermal .
By substituting the appropriate expressions for I lum (t) and θ(t) (i.e. Supplementary Equation (2) and Supplementary Equation (4), respectively) into Eq. (5) and carrying out the integration, the analytical solution for θ measured at any value of τ lum can be expressed as This analytical solution is plotted for different values of R thermal in Fig. 6(c) in the main text.

Supplementary Note 8: Confirmation of successful laser modulation up to 10 kHz
For the modulated excitation ratiometric measurements described in the main text, the output of the 980 nm laser diode was modulated using a function generator. One potential concern is that the laser output may not be successfully modulated at high frequencies.
To ensure that the modulated laser output still resembled a square wave at the highest frequencies used in our experiments, we used an avalanche photodiode (APD) to measure the 980 nm signal reflected from a silicon wafer over one excitation period for frequencies between 10 Hz and 100 kHz, as shown in Supplementary Figure 5. Above 10 kHz, the laser output begins to deviate significantly from the expected square wave shape, and we thus exclude data taken at frequencies above 10 kHz. kHz. We exclude data taken at frequencies above 10 kHz since the laser output deviates significantly from a square wave.

Supplementary Note 9: Error bars for modulated excitation ratiometric data
The error bars for the experimental data shown in Fig. 6(c) differ from the error bars for all other ratiometric data because, in this case, the ratio values are not converted to an apparent temperature rise. The final quantity that is plotted is θ measured (f exc )/θ measured (f exc = 0). Because the change in the ratio with respect to the zero-power value, i.e. ∆r = rr intrinsic (T = 296 K), varies nearly linearly with θ over the temperature range of interest, θ measured (f exc )/θ measured (f exc = 0) can be well approximated as ∆r(f exc )/∆r(f exc = 0). This  This result is not unexpected, because our objective has no built-in coverslip correction and thus cannot compensate for the refractive index contrast between glass and air 9 (a coverslip-corrected objective was not available). Upon adding a drop of water to cover the same nanoparticle that was previously imaged through the coverslip, no further loss in emission intensity was observed. We were thus able to obtain auxiliary r(I exc ) data for a particle surrounded by water. The water drop was monitored throughout the experiment to ensure that it did not evaporate. A challenge, however, is that it is difficult to determine the true excitation intensity seen by the nanoparticle. While refraction aberrations will broaden the focal spot and reduce the local excitation intensity, we measure no significant difference in the total laser power when the laser beam passes through a coverslip due the large area of our optical power meter (incidentally confirming that absorption and reflection by the coverslip are negligible). Thus, we plot the r values as a function of the power measured at the entrance to the microscope. Supplementary   Fig. 6 shows that the results for a nanoparticle in this configuration are similar regardless of whether it is surrounded by air or water, despite the fact that water has a thermal conductivity more than an order of magnitude higher than that of air. Consequently, this measurement further strengthens the conclusion associated with Fig. 4(a) and (b), which is that none of the external thermal resistors controls the apparent temperature rise. Here, we extend the results of Fig. 2 together. We measured the luminescence intensity ratio of two such double-particle clusters as a function of temperature and excitation intensity. For each double particle cluster, we also performed the same measurements on a nearby single particle.

Supplementary
Supplementary Figure 7 shows that the apparent temperature rise of the double-particle clusters is essentially indistinguishable from that of the single particles within experimental noise. From a heat transfer perspective, the fact that the apparent temperature rise of the double-particle clusters and the single particles is the same further confirms that this effect is non-thermal. The emission intensity of the double-particle clusters is exactly twice that of the single particles (see insets of Supplementary Fig. 7(a)), suggesting that two particles are clustered together in the center of the laser spot. If the two particles were further apart, yet still both within the laser spot, we would expect the emission intensity to be notably less than twice that of a single particle due to the effective decrease in excitation intensity. The self-heating estimate for a double particle cluster is thus similar to that for a single particle, with the addition of a symmetry plane at the junction of the two particles, which can be treated as an adiabatic surface 10 .
Consequently, we expect the temperature rise of a double particle cluster to be up to twice that of a single particle, yet Supplementary Fig. 7(b) shows an essentially identical apparent temperature rise in both cases. More broadly, this result suggests that the intensity-dependent photophysics we observe in this work for single particles should also be considered for ensembles.  R apparent ≈ 4.1 x 10 10 K W -1 R apparent ≈ 4.1 x 10 10 K W -1 R apparent ≈ 4.7 x 10 10 K W -1 R apparent ≈ 5.0 x 10 10 K W -1 I exc = 1.5 x 10Т W cm -2