Neuronal imaging with ultrahigh dynamic range multiphoton microscopy

Multiphoton microscopes are hampered by limited dynamic range, preventing weak sample features from being detected in the presence of strong features, or preventing the capture of unpredictable bursts in sample strength. We present a digital electronic add-on technique that vastly improves the dynamic range of a multiphoton microscope while limiting potential photodamage. The add-on provides real-time negative feedback to regulate the laser power delivered to the sample, and a log representation of the sample strength to accommodate ultrahigh dynamic range without loss of information. No microscope hardware modifications are required, making the technique readily compatible with commercial instruments. Benefits are shown in both structural and in-vivo functional mouse brain imaging applications.


Theoretical calculation of SNR and dynamic range gain
Our goal here is to evaluate the expected gain in dynamic range provided by active illumination. For this, we will evaluate the expected SNR of our system.
We begin by writing S = XP α , where S is the fluorescence signal in units of number of detected photons per sample time (i.e. unitless), X is the sample strength, P is the illumination power, and α is the multiphoton order. We also define: S 0 = set point for signal level when active illumination is engaged. S sat = signal level at which detector saturates. P 0 = laser power used for standard microscopy (w/o active illumination). P max = maximum laser power allowed when active illumination is engaged. η S = S 0 /S sat = normalized set point for active illumination. η P = P 0 /P max = normalized laser power w/o active illumination. X sat = S sat /P α 0 = maximum sample strength that can be measured w/o active illumination. B = digitization bit depth for microscope acquisition electronics. B F = digitization bit depth for FPGA input.
Our measurement of X contains noise. That is, we can write X = X + δ X, where ... corresponds to an ensemble average and δ X is the noise associated with any given measurement of X. The variance of this noise is denoted by σ 2 X = δ X 2 , leading to a SNR associated with the measurement of X given by SNR = X /σ X . This definition of SNR will be used throughout.
Let us first calculate the SNR associated with standard multiphoton microscopy (i.e. without active illumination). In this case, the digitized signal can be written asŜ = Ŝ + δŜ, with where G M is the gain associated with the microscope acquisition electronics. The first term on the r.h.s. of δŜ 2 corresponds to shot noise. To fill the bit depth of these electronics, we use G M = 2 B /S sat . Also σ 2 NS is the variance of electronic noise associated with the detection ofŜ. This includes both detector noise (in post-digitized units) and inaccuracies introduced by the digitization process itself.
Since our goal is to measure X, we should in principle also have a knowledge of the laser power P 0 . Let us assume a best case scenario where P 0 is known perfectly. We have then and to first order This may be rewritten in terms of X sat and S sat , obtaining: Standard microscopy mode: from which we readily derive the SNR associated with a measurement of X.
We proceed in a similar manner in the case of active illumination. We writeŜ = G S S andP = G P P, andL X = log 2X +C, where C is a constant introduced to ensureL X does not become negative. These quantities are digitized values processed by the FPGA board. The FPGA outputL X is then re-digitized when it is read by the microscope acquisition electronics, leading to L X = G XLX + δ N X , where δ N X takes into account inaccuracies due to digitization (i.e σ 2 N X = 0.083 for randomly distributed measurements). What is displayed on the computer screen is then L X .
To evaluate the SNR associated with this measurement, we must numerically reconstruct X. This is given by leading to For ease of calculation, we introduce here an assumption. In the case whenP is large, its measurement becomes shot-noise limited and δP/ P 1 (i.e. its contribution can be neglected). In the opposite case whenP is small, its measurement becomes electronic-noise limited, since it is detected with a simple photodiode. That is, we write δP 2 = σ 2 NP . After some algebra, we obtain Two scenarios must be considered since active illumination can operate in either feedback-active or power-limited modes.
In feedback-active mode, we have S = S 0 = η S S sat and P α = S 0 /X, leading to In power-limited mode, we have P = P max and S = P α max X = 1 η α P S sat X X sat , leading to In the event G P is properly adjusted to fill the FPGA input bit depth, then G P = 2 B F /P max and these expressions can be simplified further: Feedback-active mode: (2) Power-limited mode: We note that the transition between feedback-active and power-limited modes occurs at a threshold sample strength given by Let us now compare the maximum dynamic ranges of standard versus active illumination microscopes. With standard microscopy, the minimum possible measure of X is 1/G S P α 0 . The maximum possible measure is S sat /P α 0 . The best dynamic range, corresponding to the ratio of the two, is then G S S sat .
With active illumination, the minimum possible measure of X is obtained in power-limited mode and given by 1/G S P a max . The maximum possible measure is obtained in feedback-active mode and given by G α P η S S sat . Assuming G P is properly adjusted, the best dynamic range is then 2 αB F η S G S S sat .
The potential dynamic range gain afforded by active illumination is thus 2 αB F η S . A more conservative estimate of this dynamic range gain can be evaluated based instead on the rolloff in the SNR that occurs in feedback-active mode when X becomes large (see Fig. 1c). This rolloff occurs when uncertainties due to electronic or digitization noise in the measurement of P become larger than uncertainties due to the shot noise in S 0 , which, according to Eq. 2, occurs when . A reduced estimate for the gain in dynamic range becomes then 2 αB F η S ξ . This reduced estimate is likely overly conservative since the SNR at X r still remains high.
At this point, it is instructive to attach numbers to our calculations. Parameters needed to establish SNR as a function of normalized sample strength X /X sat are η S , η P , G S , G P , and S sat . As an example (c.f. Fig. 1c), we use η S = 0.5 and η P = 1, where the latter ensures that photodamage cannot increase. We also assume the electronic gains are adjusted such that G S = 2 B F /S sat and G P = 2 B F /P max (bearing in mind that we could do slightly better with G S = 2 B F /S 0 ). Finally, we need to calculate S sat in units of number of detected photons per sample time. A typical PMT sensitivity is 10 5 A/W and a typical PMT output saturation level is 100µA. Hence PMT saturation typically occurs at roughly 1nW, or 2.5 × 10 9 photons/s. Assuming a sampling time of 1µs and a PMT quantum efficiency of order 10%, we arrive at S sat ≈ 250.

FPGA-based PID design
Our feedback system is designed based on the Red Pitaya development kit. Red Pitaya is a FPGA development kit with a dual core ARM Cortex chip on board and features two 125 MS/s analog inputs and two 125 MS/s analog output ports. The hardware parameters and functional layout are summarized below. As is illustrated in Fig. S1, our FPGA-based PID controller is composed of five functional blocks.  Figure S1: Overall structure of our PID controller. Functional blocks are shaded in gray. s(k) is the discretized PMT signal; s 0 is the set point; k P and k I are feedback coefficients; u is the controller output signal, which controls the transmission of the EOM; p is the laser power, measured by the photodiode; x is the sample strength.
1. A PID controller that computes the desired controller output.
2. An AXI-Stream interface for the user to set the feedback parameters.
3. An oscilloscope to visualize the input signal.

A sample strength calculator.
5. An input/output interface.
Among them, the advanced extensible interface (AXI), the oscilloscope and the I/O interface is part of the Red Pitaya open source project. Detailed code and description can be found on github: https://github.com/RedPitaya/RedPitaya/tree/master/fpga The AXI-Stream interface provides real-time communication between Red Pitaya's Linux system and FPGA. We use it to adjust the feedback parameters and visualize the input signal. The input/output interface performs dual channel A/D and dual channel D/A conversion at 125 MS/s. In our case, signal level is digitized by input channel 1 while the illumination power is digitized by channel 2. The oscilloscope modules stores data from the high speed A/D converter in a triggered buffer and transfers the data through AXI-Stream to Red Pitaya's DRAM for on-screen visualization.

4/9 PID controller
Our PID algorithm is a widely-used negative feedback control algorithm that begins with a calculation of a time-varying error signal, e(t), which is the difference between the values of the set point and current signal level. The error signal is sent to three contributing control parameters: a control directly proportional to the current error, a control based on the integrated past error signal, and a control based on the derivative of the error signal. The continuous PID control algorithm in standard form can be written as (12) In our case, we did not use the derivative component of the algorithm and set this to zero. As shown in Figure S2, we first subtract the input signal from the set point to obtain the error signal. We then multiply the error signal by a proportional gain k p to obtain the proportional control term. At the same time, the error is multiplied by an integral gain k I and added to a integral register. Finally, the proportional and integral terms are summed and the result is sent to an output actuator. The implementation of the algorithm contains two multipliers, two adders and one subtractor. To improve calculation speed, Xilinx R IP cores are used.

PID Controller
Register -e(k) Figure S2: Work flow of the PID controller.
To alleviate the trade-off between feedback response time and ringing, we regulate the feedback gains in real time. In brief, when the imaging system scans through a dark area, the feedback gain is increased to reduce its response time; when the system scans through a bright area, the feedback gain is reduced to mitigate ringing.

Sample strength calculator
The use of a logarithm scale to encode sample strength was introduced in Chu et al (4). As shown in Figure S3, three steps are involved. First, the signal level is measured by the PMT and digitized as S and the illumination power is measured by the photodetector and digitized as P. Both S and P are encoded as 13-bit signals, meaning their logarithms can be computed using a same predefined and scaled lookup table that occupies less than 7KB of on-board ROM. Finally, a subtraction operation is performed between the logs of S and P (the latter being multiplied by 2 in the case of two-photon microscopy by bit shifting). A constant is added to the final result (here 2 13 ) before sending it to the microscope control system to prevent it from going negative.

FPGA Implementation
FPGA programming is performed through the Vivado R Design Suite 2015.1 by Xilinx. All source code is written in Verilog. The target FPGA device used in this research is the Xilinx Zynq 7010 SoC.

Device Resource Utilization
The device resource utilization is shown in Table 2.

Timing Analysis
There are several devices in our feedback loop, including a PMT, its associated transimpedance amplifier, and the EOM driver. Altogether, these devices were found to introduce a latency in the feedback of ≈ 2µs. In contrast, the delays introduced by the FGPA operations equipped with a 125MHz oscillator as the master clock are: 1. A/D conversion. For the LTC2145 A/D converter, ≈ 7ns.

D/A conversion. For Red Pitaya's on-board D/A converter, ≈ 30ns.
The total latency of ≈ 50ns introduced by the FPGA is thus negligible compared to that introduced by the devices. In our case, the dominant speed limitation of our feedback comes from our EOM driver.
In our experiments we used a pixel sampling rate up to 1MHz ( Figure S4 and Supplementary Video 9). Such a high pixel rate pushed our system to its very limit since barely a single control cycle is allowed per pixel. With our gain setting, the PID control can bring the PMT signal into a measurable range in the first cycle and thus avoid saturation. However, when only one cycle is allowed per pixel problems inevitably occur in PID feedback such as overshooting on upswings and undershooting on downswings. These problems are highlighted in Figs. S4 and S5, which illustrate active illumination at 1.2µs and 4µs pixel times. The sample strength log 2 X, the raw PMT signal S, and illumination power P are simultaneously acquired with three channels of the microscope acquisition electronics. Ideally, the raw PMT signal should plateau to a constant value whenever the system toggles to feedback-active mode. Manifestly, even at 4µs pixel time, some overshooting occurs before the system is able to stabilize. Nevertheless, the calculated sample strength X remains unaffected because transient errors in S and P largely cancel one another. 6/9 Figure S4: Images of YFP-labeled neurons in neocortical slice at 1.2µs pixel time (833kHz pixel rate). S is the raw PMT data while active illumination is on. P is the illumination power measured by the photodiode. log 2 X is the resultant log-encoded sample strength sent to the microscope acquisition electronics. X is the linear-encoded sample strength calculated post-hoc. Note the bright edges apparent in S are caused by overshoot in the feedback at 833kHz pixel rate. These are canceled after computation of log 2 X (or X). Reduced photobleaching with active illumination Because AI effectively sets a limit on the fluorescence signal level when it is active, we expect it to also set a limit on the rate of photobleaching. To demonstrate this, we imaged a dendrite from a YFP-labeled cortical pyramidal cell. The total laser power incident on the slice was about 20 mW. We performed a sequence of 20 volumetric scans over the same region, where AI was ON for the first 10 scans and OFF for the remaining 10 scans. When AI was ON, no saturation occurred. When AI was off, occasional isolated spots in the dendrite caused saturation. We selected of region of interest (circumscribed in red) that avoided these spots. The total fluorescence signal produced from within this region is plotted as a function of time (or scan number). Manifestly, when AI is on, the photobleaching occurs at a slower rate than when AI is off, confirming the unsurprising result that the photobleaching rate becomes reduced with the reduced excitation power occasioned by AI. Figure S6: A YFP-labeled dendrite repeatedly imaged with (a) AI on and (b) AI off. The darker regions in (c) correspond to regions where AI became active. The total fluorescence power recorded from a region of interest (circumscribed in red) is plotted as a function of scan number. The plots are normalized to their respective starting levels. When AI is ON, the photobleaching rate is slower than when AI is OFF.