T ENSORIAL TOMOGRAPHIC F OURIER P TYCHOGRAPHY WITH APPLICATIONS TO MUSCLE TISSUE IMAGING

,


Introduction
Quantitative phase imaging (QPI) is a well-known label-free microscopy approach that can detect phase delay introduced by semi-transparent cells and tissue [1].Due in part to its ability to provide quantitative information about primarily transparent biological specimens with low phototoxicity, QPI has become an invaluable tool in scientific and clinical studies, including for monitoring neuronal firing [2] and cancer cell line detection [3], to name a few.Besides scalar phase contrast, transparent specimens also exhibit alternative and important endogenous optical contrast mechanisms, including anisotropic properties such as material birefringence and orientation.Indeed, the orientation of molecular arrangements in lipid membranes can now be quantitatively monitored with polarization-sensitive microscopes for studying multi-organelle interactive activities [4,5].There are naturally a variety of polarizationsensitive microscope arrangements, including early analog designs [6,7], differential interference contrast (DIC) methods [8], and contemporary digital approaches that reconstruct quantitative specimen retardance and orientation [9,10,11].These methods have been applied to study small model organisms [12,13] and to assist clinical diagnosis [14,15,16,17].
Recently, there has been increased interest in measuring polarization-sensitive phase information from specimens.In general, polarization-sensitive quantitative phase imaging (PS-QPI) methods can be divided into two categories: those that rely on interferometric detection, and those that utilize computational phase retrieval methods.Interferometric methods (e.g., off-axis holography) can create polarization-sensitive phase images with as few as a single measurement [18,19,20,21], and can be extended to 3D with diffraction tomography approaches [22,23,24,25].While often impressive, these methods usually require complex arrangements of coherent laser illumination and careful system alignment, which sets additional design requirements for use in clinical applications.Computational phase retrieval methods, on the other hand, rely on image reconstruction algorithms that convert multiple intensity measurements into phase-sensitive outputs [26,27,28,29,30,31], and can be implemented with less expensive hardware.Due to their simple instrumentation, computational polarization microscopes have been increasingly applied to image biological samples, for instance, to study white matter tracts within whole brain slices [30] and to diagnose malaria from blood smears [28].Furthermore, these methods have been extended to image 3D samples, such as axons and cardiac tissue with axial scanning [32,33].Thick samples such as organoids and tissue slices also naturally have intriguing 3D structures that can include anisotropic material.Creating high-resolution volumetric representations of the polarization properties of these samples is essential to studying biology and pathology.While tomographic imaging methods such as confocal-based approaches have been developed in the past to image muscle tissue and neural organoids, for instance [34,35], there remain relatively few microscopic techniques to jointly capture quantitative phase and anisotropy across a large three-dimensional volume at high resolution.One very recent study attempted to create 3D anisotropy maps using off-axis LED illumination but did not provide tomographic permittivity matrix reconstructions [36].
Here, we propose a non-scanning polarization-sensitive tomography method, termed Tensorial tomographic Fourier Ptychography (T 2 oFu), to create quantitative volumetric permittivity matrix reconstructions without any moving parts.Our method is an extension of recently developed intensity optical diffraction tomography principles [37,38,39,40,41,42] and has the potential to be extended to a video-rate system in the future [37].We image a variety of calibration targets, as well as 3D tissue and muscle fibers samples, through which we demonstrate the ability to resolve the fine filament structure and are consistent with one of the gold standards, i.e., second harmonic images, found in the literature [34].Additionally, we demonstrate reconstructions of a heart tissue sample that carries important information for detecting cardiac amyloidosis.

Experimental setup
The imaging setup of T 2 oFu is illustrated in Fig. 1(A).The illumination consists of an array of 25 × 25 addressable LEDs (µInventions Inc; Portugal) and a left-circular polarizer (CP42HE; Edmund Co., USA).521 nm wavelength light from the LED array is circularly polarized and illuminates the 3D sample from different angles.A microcontroller (ARM Cortex-M3) and a voltage level shifter (SN74AHCT) are used to turn on each small LED pixel (WS2812b-2020) sequentially.The optical field then passes through an optical system consisting of an infinity-corrected objective (0.25 NA or 0.4NA, Plan Achromatic; Olympus, JP) and a tube lens (Achromatic Doublets, 180mm focal length; Thorlabs, Inc, USA).The largest illumination NA (N A illum ) is chosen to match the native NA of the objective (N A obj ).Hence, the theoretical lateral and axial resolutions are δ x,y = λ /2(NA obj +N A illum ) and respectively [43].For the 0.25 NA system, δ x,y = 5.2 µm and δ z = 8.2 µm.For the 0.4 NA system, δ x,y = 3.2 µm and δ z = 3.1 µm.The intensity image is captured with a polarization-sensitive CMOS camera (BFS-U3-51S5P; Teledyne FLIR LLC; USA).The polarization camera analyzes the light orientated at four different polarizations (0°, 45°, 90°, and 135°), achieved by placing 2×2 wire-grid polarizers between the pixel array and the microlens array, as illustrated in the m-dashed box in Fig. 1(A).Exemplary intensity images under illumination from different angles are shown in the box in Fig. 1(B).

Notation
Here we introduce the notation we use in the rest of this article.First, we use • and • symbols to denote vectors and matrices for variables, respectively.By default, all the vectors in this article are column vectors.Second, we utilize the Mathematical Script font enclosed in curly braces to denote spatial operators.As an example, one frequently used operator is F{•}, which indicates the Fourier transform in space.Next, bold letters in lowercase represents support vectors in either frequency (u) or space (r).Finally, we utilize the • symbol to denote the frequency-domain counterparts of variables previously defined in the space domain, such as Ĩ(u) = F{I(r)}.

Vectorial light propagation
The optical properties of a 3D sample can be described by its 3 × 3 permittivity matrix [44] where r = (x, y, z) is the voxel position in space.In the scope of this work, we adopt the first Born approximation [32,23], and the scattered vectorial electric field E s and illumination E 0 is related as where V (r ) = ¯ − ¯ 0 is the sample scattering potential tensor with ¯ 0 as the permitivity tensor of the background medium.Ḡ(r) is the dyadic Green's tensor [45].In addition, we use illuminations with relatively small incident angles, which have weak polarization along the z-axis.Therefore, we mainly consider the transverse polarization of the electrical field, E ⊥ = (E x , E y ) T , which will be discussed in more detail in the next subsection.The intensities of the vectorial electric field E ⊥ at different orientations are then analyzed and recorded by a polarization-sensitive optical imaging system, which we model as a 2 × 2 pupil matrix at each spatial frequency u [26].Therefore, the measured intensity under plane wave illumination with angle u analyzed by a polarizer with Jones vector a is For a linear polarization analyzer oriented at α, a = [cos α, sin α] T [46,26].The intuition behind this is that light can only oscillate in the same direction as the linear polarizer [26].

Forward model and inverse problem
For our initial demonstration of T 2 oFu, we modify the model discussed in Section 2.2.2 by making certain approximations.The objective of this is to craft a meaningful forward model that accurately describes our experimental measurements, whose inverse problem is less ill-posed.First, following Saba et al. and other closely relevant works in literature [23,27,28,26], we begin with a paraxial approximation, which assumes weak polarization along the optical axis of the illumination and negligible interaction between traverse and axial polarization from the sample.With this, we can simplify the 3 × 3 permittivity matrix to a 2 × 2 matrix, denoted as While this approximation may not always be accurate for every anisotropic sample when illuminated at high angles, it is accurate under up to a 25°oblique illumination, based on a study using finite element analysis [23].Additionally, we assume that the background media is isotropic and uniform (i.e., not spatially dependent) with a diagonal permittivity tensor ¯ 0 = 0 I.This simplifies Green's tensor to a diagonal matrix with the same component for each polarization [23], G(r, r ) = G(r), where r = |r − r |, is the scalar Green's function that has a Weyl expansion [44,47] with wavenumber vectors in lateral (u = (k x , k y )) and axial (η = k 2 0 − |u| 2 ) directions.k 0 is the wavenumber of the isotropic background medium.This expansion is generally easier to work with when we prefer to represent the object in both lateral frequency and axial space domain.
In addition, we make the assumption that the sample being imaged is homogeneous, which is usually assumed for many types of crystals and biological samples [48,10].The permittivity matrix is then symmetric and can be decomposed into [48]: The variables o and e represent the permittivity values along the ordinary and extraordinary axes, respectively.The parameter θ is the angle between the principle axis and the extraordinary axis (also known as the slow axis).Following the convention introduced by an early literature [49], we rename elements in Eq.4 as Since we use left circularly polarized illumination, the transverse scattering potential becomes where Finally, we refer to Ling et.al. [50] for two further approximations of the illumination and scattering processes.Our first assumption is that the illumination from each LED at the sample plane is a plane wave, which is commonly used in Fourier ptychographic computational microscopy [51].The second approximation involves utilizing a weak object assumption that ignores the second-order scattering term [37,52,53].Further, we disregard pupil aberration and model it as a low-pass filter P (u) with a cutoff frequency based on the numerical aperture of the objective lens [37].Jones matrix has been demonstrated feasible in previous literature for the correction of anisotropic aberration by jointly reconstructing the pupil [26,31], which is planned for future research.Based on these approximations, the forward model can be expressed as [37] Ĩl,m (u, where Ĩl,m (u, z) and Ĩl,m 0 (u, z) are 2D Fourier transform of the measurement and DC term from m th LED illumination analyzed by l th polarizer, respectively. and are the diffractive transfer functions in frequency and space for the real and imaginary part of the scattering potential under m th LED illumination with shape S m (u ) [52,54].Since each LED has a very small die area (< 170 µm in diameter), we represent it as a delta function.l ∈ {0°, 45°, 90°, 135°}.Ṽ l Re (u, z) and Ṽ l Im (u, z) are 2D Fourier transforms of V l Re (r) and V l Im (r) along lateral directions at depth z, respectively, and are related to Eq. 10 via V 90°R e (r) + jV 90°I m (r) = V 3 (r) + jV 2 (r) which are the scattering potential components corresponding to each analyzer angle.We want to point out that V 1,2,3 (r) are complex variables that may have imaginary parts; hence, the above equation does not imply ṼRe (r) = V 1 (r), for example.
For concise expression, we define a new variable v ∈ C N ×M ×T ×3 representing all the potentials V 1 (r), V 2 (r), and V 3 (r), where N, M, T are width, height, and depth of the 3D sample, respectively.Further, we introduce the operator A l,m (•) as the forward model for the sample illuminated by the m th LED and analyzed by the l th linear polarizer.To reconstruct the permittivity matrix, we formulate the inverse problem as with the loss function tv(•) is the isotropic total variation operator.γ is a regularization parameter empirically set to be 1 × 10 −6 for all experiments.The forward model is implemented in Pytorch and the loss function is optimized using a stochastic gradient descent method with Nesterov momentum acceleration [55,56].
Subsequently, we extract polarization properties of interest, such as orientation and birefringence from T 2 oFu reconstructions 1,2,3 .The permittivity of ordinary and extraordinary axes can be computed as and the refractive index along the ordinary and extraordinary axes n o,e = √ o,e can be further derived, along with the averaged refractive index (n = n o + n e ) and birefringence (∆n = n e − n o ).Moreover, the orientation to the slow axis is computed as

Experimental Results
To validate the performance of the proposed method, we first show results from a variety of calibration targets.Due to the lack of commercially available high-resolution 3D polarization phantom, we follow previous works to validate different aspects separately [32,28,57,26].We then show reconstructions of a single fixed muscle fiber.Finally, we show reconstructions of a tissue section sample from a human heart biopsy that is predictive for lethal cardiac amyloidosis.All the images presented here are captured and then reconstructed using a 20×, 0.4 NA system, except for the cardiac tissue, which is recorded and processed with a 10×, 0.25 NA system.

Tomographic reconstruction
In this subsection, we demonstrate tomographic reconstructions of averaged refractive index n and birefringence ∆n from isotropic and anisotropy calibration samples.These results are acquired with a 20× objective (0.4NA, Olympus Corporation, Japan.)imaging system.Figure 2 (A) shows the image captured with center LED illumination when the imaging system is focused at the middle of a polystyrene microsphere sample.The sample is made of two layers of 3µm-diameter microspheres immersed in n = 1.575 oil.(B) shows images taken with center LED illumination when focused at different depths.The arrows in (B) highlight the microspheres that are in focus, suggesting two layers of microspheres are placed at −7µm and 5µm planes.Note that since the refractive index of polystyrene (n = 1.60 @ 520nm) is very close to the background medium, the contrast of in-focus microspheres is very low.(C)-(D) shows the reconstructed refractive index and birefringence at two different depths.The accuracy of the tomographic depth reconstruction is validated with images displayed in (B).Since polystyrene is an isotropic material, the reconstruction shows no anisotropy properties except on the edge, which agrees with the well-recognized edge birefringence phenomena [58] reported in prior literature [32,26].(E) shows the cross-sections of the reconstructed reflective index and birefringence.(F) plots the profile of the reconstructed birefringence and reflective index averaged across 10 microspheres.
In addition, we show reconstructions of an anisotropic potato starch sample in Fig. 3.The sample consists of two potato starch grains immersed in n = 1.515 oil at different depths.Fig. 3(A) depicts captured intensity images illuminated with the center LED at four different polarizations.We can see the spiral patterns imaged with different polarization differs from each other noticeably.Similarly, Fig. 3(B) shows images when the system is focused at different depths of the sample.The black arrows point at grains that are in focus, suggesting the two potato starch grains are suspended at −9 µm and 6 µm.Fig. 3(C) shows the reconstructed orientation and birefringence.To best visualize the results, we follow the convention [21,28,26] to display this multidimensional data using an HSV colormap, where saturation is set to one, value is associated with birefringence, and orientation is coded in hue.The reconstructed structures agree with starch grain reconstructions reported in previous holography-based literatures [23,25].

Orientation measurement
To verify the accuracy of the reconstructed orientation, we test our method on a sample made with monosodium urate (MSU).MSU are needle-shaped crystals precipitated from uric acid that could trigger robust inflammation such as acute arthritis and other immune activations that cause severe pain in patients [59].Figure 4(A) displays the reconstructed birefringence and orientation at two slightly different depths.These results were also obtained with a 20× objective imaging system and suggest that the reconstructed orientation values of line-shaped MSU crystals follow the structural direction (labeled next to each MSU crystal).This is in agreement with reconstruction results from LC-PolScope-based methods found in the literature [28,57].Additionally, Figure 4(B) presents a 3D rendering of the reconstructed birefringence.

Muscle fiber assessment
High contrast and high-resolution structural imaging of intrinsic signals in muscle skeletal fibers is an important task for the rapid detection of changes in myofibrillar organization that can lead to skeletal myopathies.Currently, 3D muscle tissue is typically imaged by complex and expensive systems, such as second harmonic generation (SHG) microscopy.SHG exploits the contrast of polarization properties from the non-linear susceptibility in molecules like myosin, by using high-power, ultra-short pulsed lasers in a point scanning configuration [34].Here we show 3D T 2 oFu reconstruction results of an isolated healthy muscle fiber using our inexpensive, LED-based and scanning-free system.Figure 5(A) shows an image of a muscle fiber captured with center LED illumination.region of the image in x-y resembles reconstructions reported in [34].Figure 5(C) shows an image of the same muscle fiber from a different field of view, where a non-muscle fiber with a ninety-degree bend is placed below the muscle fiber.pattern of a healthy muscle fiber.The zoom-in region highlights the muscle grains with consistent orientations, in agreement with the results reported by Both et al. [34].The reconstruction shows a change in orientation at the bend of the non-muscle fiber (D1), while the orientation of the muscle fiber remains constant.Furthermore, the spatial resolution of the reconstruction is sufficient to clearly resolve the muscle filaments (D3).

Imaging cardiac amyloidosis
Finally, we apply our method to image a heart tissue sample that has cardiac amyloidosis.Cardiac amyloidosis is a lethal disease that affects more than 12,000 patients in the US alone with < 5% 10-year survival rate [60].In current practice, the biopsied tissue is first frozen and thinly sliced, then stained with a congo red-colored dye and inspected under a cross-polarized microscope.Figure 6(A) shows the brightfield image.(B) shows the same region, imaged with a cross-polarized color microscope.The vibrant apple green color suggests mostly likely amyloid protein has built up inside the tissue sample.(C)-(D) depict the reconstructed refractive index and birefringence, while zoom-ins of the boxed region at different depths are shown in (E)-(F).Since the sample is thinly sliced, we do not observe noticeable structural changes in different layers.However, the structure of the birefringence reconstruction is correlated with the color-stained cross-polarized image, which could potentially be useful for rapid on-site inspections in the future.

Discussion and Conclusion
In this article, we introduce T 2 oFu, a new non-scanning microscopy method that reconstructs volumetric permittivity metrics of samples based on computational illumination strategies to record polarized measurements and retrieve phase information.Using relatively low-NA objectives, we demonstrate that T 2 oFu can provide polarization-sensitive tomographic reconstruction for various calibration samples and biological specimens that are potentially useful for future scientific and clinical studies.

� 𝑛𝑛 Δ𝑛𝑛
Cross-polarized Brightfield To ensure successful clinical translation, there are a few improvements that can be made in the next step.First, to increase the frame rate of T 2 oFu, high flash-rate LEDs need to be deployed [61].Novel sensitive camera sensors [62,63,64] will also play important roles in maintaining a good signal-to-noise ratio.In terms of the reconstruction algorithm, as a first demonstration, we use a gradient-based method implemented with auto-differentiation to optimize the loss function in Eq. 16 [55].As the data fidelity term has a closed-form solution, variable-splitting methods can be very effective and allow incorporating advanced regularization that does not have implicit forms [65,66].Moreover, due to spectral scanning, so far, T 2 oFu has a lower frame rate compared to most holographic methods in the literature.To improve T 2 oFu's throughput, in addition to using higher-speed LEDs and a more sensitive camera [37], adopting a multi-aperture approach appears to be another promising direction [67,68,69,70,71].Finally, concerning modeling, to relieve the ill-posedness of the problem, we have approximated the permittivity tensor with its lateral components in this work.Albeit widely used [23,27,28,26,72,73,74], this simplification disregards out-of-plane anisotropy.Recent works have shown under uni-axial approximation, tomography of 3D polarization orientation can be retrieved [32].Further investigation in this direction to extend the current method to extract valuable out-of-plane information is planned [32,75].

Figure 2 :
Figure 2: Reconstruction results for polystyrene microspheres.(A-B) Intensity images of the sample illuminated with center LED.(A) The data is captured when focused at z=0 µm.(B) Images were captured when focused at -5 µm and 7 µm by mechanically moving the sample.These images serve as a reference to be compared with the reconstruction.(C)-(E) plots the reconstructed refractive index and birefringence.All columns in (C) and (D) share the same scale bar, all the reflective index reconstructions (C), as well as all birefringence reconstructions (D) and all cross-sections (E) share the same colorbar.(E) displays cross-sections of the reconstruction in places color labeled in (A).(F) plots the profile of reconstruction averaged over 10 microspheres.

Figure 5 (Figure 3 :
Figure 3: Reconstruction results of potato starch grain.(A) Polarized intensity measurements from the center LED.(B) Images taken when focused at -9 µm and 7 µm serve as a reference for validating the reconstruction.(C) Tomographic reconstruction of the sample.The orientation is coded in color, while the birefringence is represented as brightness.All the images in (C) share the same colorbar.

igure 4 :
Reconstructions of monosodium urate (MSU) crystals.(A) shows the reconstructed polarization orientation of the crystals at two different depths, labeled with their respective structural directions.(B) renders the birefringence reconstruction of the sample in 3D.

Figure 5 :
Figure 5: Reconstructions of a muscle fiber.(A) shows the image of a muscle fiber with the center LED illumination.The imaging system is focused in the middle of the muscle fiber.(B) plots the reconstructed birefringence.The zoomin region resembles reconstruction results reported in[34] using a second harmonic generation imaging approach.(c) shows the image of the same muscle fiber at a different region, where a non-muscle fiber with a ninety-degree bend is placed below the muscle fiber (see -10µm in D).The imaging system is focused between this and the muscle fiber.(D) shows the reconstructed orientation at different depths, with a zoom-in showing the fine sarcomere structure of muscle tissue.

Figure 6 :
Figure 6: Images of a heart tissue sample with cardiac amyloidosis.(A) displays a brightfield image.(B) shows a cross-polarized image taken with a color microscope.(C)-(D) present reconstructed refractive index and birefringence, along with zoom-ins of the boxed region at different depths depicted in (E)-(F).(C) and (E) share the same color bar, while (D) and (F) share another common color bar.