Phase-contrast X-ray tomography resolves the terminal bronchioles in free-breathing mice

Phase-contrast X-ray lung imaging has broken new ground in preclinical respiratory research by improving contrast at air/tissue interfaces. To minimize blur from respiratory motion, intubation and mechanical ventilation is commonly employed for end-inspiration gated imaging at synchrotrons and in the laboratory. Inevitably, the prospect of ventilation induced lung injury (VILI) renders mechanical ventilation a confounding factor in respiratory studies of animal models. Here we demonstrate proof-of-principle 3D imaging of the tracheobronchial tree in free-breathing mice without mechanical ventilation at radiation levels compatible with longitudinal studies. We use a prospective gating approach for end-expiration propagation-based phase-contrast X-ray imaging where the natural breathing of the mouse dictates the acquisition flow. We achieve intrapulmonary spatial resolution in the 30-μm-range, sufficient for resolving terminal bronchioles in the 60-μm-range distinguished from the surrounding lung parenchyma. These results should enable non-invasive longitudinal studies of native state murine airways for translational lung disease research in the laboratory. Mechanical ventilation of living animals is routinely used to achieve high-resolution pulmonary imaging, but this can damage the subject. Here, an alternative, free-breathing method enables X-ray tomography with 30 μm resolution.

I maging is an indispensable tool for functional and structural studies of the lungs as well as understanding the progression of respiratory diseases. Here, X-ray imaging and computed tomography (CT) are the key modalities for high-spatialresolution imaging. In recent years, phase-contrast X-ray imaging has broken new ground in preclinical respiratory research as it enables improved contrast at air/tissue interfaces due to the large differences in electron density. While the grating-based version 1 has attracted attention for its dark-field contrast being sensitive to the lung parenchyma microstructure (as demonstrated on animal models [2][3][4][5] and in a recent clinical study 6 ), the propagation-based technique 7,8 provides a conceptually simpler yet more dose-efficient approach to achieve phase-contrast 9 . Propagation-based phase-contrast X-ray imaging has been implemented using both synchrotron and laboratory X-ray sources and is routinely used for preclinical lung research (e.g., mice and rats [10][11][12] , rabbits 13 , and pigs 14 ). Laboratory translation (i.e., small-footprint imaging arrangements 15 ) is especially important for adoption of these techniques within the broader preclinical research community.
The major challenge for in vivo tomographic imaging of the lungs is cardio-respiratory motion. High-resolution imaging of animal models typically require acquisitions over many breathing cycles. To minimize the associated image blur, intubation and mechanical ventilation is often employed to synchronize image acquisition with the cardiac and/or respiratory cycle. However, the relations between ventilation parameters (e.g., tidal volume, positive end-expiratory pressure) and their effect on lung function are complex, rendering the prospect of ventilation induced lung injury (VILI) 16 a confounding factor when the aim is to study the lungs in their native state. Notably, mechanical ventilation in mouse models has been associated with inflammatory responses 17 , pulmonary coagulopathy 18 , increased airway volumes 19 and alveolar collapse 20 . Furthermore, intubation and forced breathing adds an invasive layer to the imaging workflow which increases the risk for complications and harm to the animal. While prospective gating on free-breathing mice as an alternative is not a new concept [21][22][23][24] , it has yet to provide the spatial-resolution matching implementations using forced breathing, i.e., resolving the smallest airways.
In the present work, we show that propagation-based phasecontrast X-ray tomography on free-breathing mice enables noninvasive and label-free 3D visualization of the tracheobronchial tree, achieving intrapulmonary spatial resolution down to the 30-μm-range with our small-footprint laboratory arrangement. Terminal bronchioles in the 60-μm-range can be distinguished in 3D at radiation doses (300-450 mGy) and scan times (12-20 min) allowing for long-term respiratory studies of the same animal (i.e., longitudinal imaging). We quantify the spatial resolution using both a micro-CT phantom and directly on in vivo tomographic data using Fourier-domain methods. Lastly, we discuss the potential impact of imaging free from mechanical ventilation for translational studies of respiratory disease, as well as developments towards alveolar-level imaging in the laboratory.

Results
Our small-animal imaging arrangement is shown in Fig. 1a. The liquid-metal-jet micro-focus X-ray source provides the small emission spot necessary for achieving propagation-based phasecontrast. A pneumatic shutter in front of the X-ray source is used to control the X-ray exposure. Animals are positioned vertically in a custom-made 3D printed support placed on top of a motorized rotation stage and kept under anesthesia during tomographic acquisition. By introducing the distances R 1 and R 2 , X-ray projections of the mice are magnified and captured on a large field-of-view mammography detector. Free-space propagation of transmitted and refracted X-rays through the animals results in phase-contrast fringes on the detector plane (cf. Fig. 1b).
We use a pneumatic pillow sensor to monitor the respiratory cycle for prospectively gated image acquisition (cf. Fig. 1c). Each breath triggers a transistor-transistor logic (TTL) pulse, followed by an adjustable delay (typically 300-500 ms) to ensure fully relaxed respiratory motion prior to image acquisition. To minimize the radiation dose, the camera exposure extends before and after the shutter opening window, which then effectively determines the X-ray exposure time (typically 100 or 150 ms). We opted out from synchronizing with the cardiac movement as well (i.e., cardio-respiratory motion) as this would impose much shorter exposure times (<10 ms) and thus require higher power X-ray sources (e.g., synchrotrons 12 ) to achieve sufficient photon statistics in each projection image. We note that our shortest exposure imaging already operates at the lower-end of the detector dynamic range (cf. Fig. 1d). Following X-ray exposure, the image is read out from the detector and the animal is rotated for a total of 720 projections equally spaced over 360 ∘ resulting in a tomographic scan. The total scan time is determined by the average breathing-rate of the animal. Our acquisition works up tõ 60 Hz, resulting in tomographic scans as quick as 12 min (i.e., 60 projections/min). The resulting radiation dose is 300 or 450 mGy per tomography depending on the choice of 100 or 150 ms X-ray exposure times per projection. We note that these dose levels are compatible for repeated exposure and therefore longitudinal studies of the same animal 25,26 .
To evaluate if our gating implementation reduces respiratory motion blur, we acquired a set of posteroanterior chest projection images of a mouse with and without gating, shown in Fig. 2. Compared to sequential acquisition without gating (cf. Fig. 2a) our prospective gating approach successfully removes the bulk respiratory motion (cf. Fig. 2b). We note that as the acquisition is triggered on the macroscopic breathing movement, microscopic intra-pulmonary motion both within and between X-ray exposures is likely still present. In addition, respiratory gating on freebreathing mice relies on the lung relaxing to the same position following each individual breath, which can further explain the residual movement observed in Fig. 2b. By tracking a single rib in between consecutive projection images, we observe that the 2D displacement is in the order of~20-30 μm or~2-4 pixels in the projection domain (cf. Supplementary Movie 1). As the latter is a 2D projection, we expect the real 3D displacement of the ribs to be larger. No differences in motion-blur reduction could be observed when adjusting the exposure delay between 200 and 600 ms at 45-60 Hz average respiration rates.
Next, we determined the spatial resolution of our imaging arrangement using a micro-CT phantom. The results are shown in Fig. 3. In the projection domain ( Fig. 3a-c) we can resolve 50 lp/mm patterns in the horizontal direction (cf. Fig. 3b), corresponding to full-period features of 20 μm (d min ) and a pointspread function (PSF) with a full width at half maximum (FWHM) in the same range. In the vertical direction, however, the corresponding frequency pattern is not as clearly resolved, explained by the elliptical X-ray emission spot (~8.5 and 18 μm FWHM in the horizontal and vertical sample plane direction, respectively) combined with the detector PSF (~13 μm FWHM in the sample plane). The larger vertical emission spot enables increased source power, important for short exposure times while preserving high resolution in the horizontal direction, which is critical for the axial tomographic resolution. While Fourier Ring Correlation (FRC) 27 can be used for an observer-free assessment of the spatial resolution (cf. Fig. 3c), the results in the projection domain exhibit artifacts both due to the elliptical X-ray spot (two intersections with the threshold-curve) and the fact that the expected resolution is too close to the Nyquist sampling (pixel size = 8.25 μm). Nevertheless, the results show that the projection imaging is detector-limited and source-limited in the horizontal and vertical direction respectively.
Similarly, we determined the axial tomographic spatial resolution ( Fig. 3d-f). Here, the modulation of the 33.3 lp/mm patterns could be resolved while the 50 lp/mm patterns were not (cf. Fig. 3e). This is in agreement with FRC (cf. Fig. 3f) predicting ã 30 μm resolution limit (d min ). In contrast to the projection domain, the monotonically decreasing FRC curve near the threshold intersection indicates that the tomographic axial resolution is isotropic. The main limiting factors of our tomographic resolution are the limited angular sampling (due to radiation dose constraints) and the phase-retrieval processing which improves contrast by low-pass filtering.
Resolving terminal bronchioles in free-breathing mice. The spatial resolution of our arrangement proves to be sufficient for observing the smallest airways in free-breathing mice. Figure 4 shows an example axial lung slice from a representative 450 mGy tomographic scan (whole lung-volume stack shown in Supplementary Movie 2). Within the lungs, areas of darker intensity represent air-filled compartments. At closer inspection we observe airways down to the 60-μm-range (cf. Fig. 4c) clearly distinguished from the surrounding lung parenchyma. Based on diameters alone, these should be terminal bronchioles (50-100 μm 28-30 ) although expected diameters are influenced both by strain and inflation-state. We note that compared to images of murine lungs acquired at end-inspiration (i.e., with mechanical ventilation 11,12 ) the lung parenchyma here appears overall brighter (i.e., higher average density) and with less contrast at, e.g., the pleura. This is expected as inspiration mainly expands alveoli rather than bronchioles 31 and results in a larger fraction of air in the parenchyma, reducing the average density. Somewhat unexpectedly, however, is that the denser parenchyma in the end-expiration state in fact facilitates studies of the tracheobronchial tree (compared to the end-inspiration state) by increasing the parenchyma/bronchiole contrast. Diameters of murine alveoli reported in literature should be in a similar range as the terminal bronchioles (40-90 μm 29,32,33 ), but being tightly packed as well as more sensitive to motion blur from residual inflation dynamics should explain why our 30-μm-range spatial resolution in the lungs (cf. air/wall interface shown in Fig. 4b) does not allow for separation of individual alveoli. A fewfold higher spatial resolution (sub-10-μm) possibly combined with shorter exposure times would be required for discerning alveoli in the parenchyma regardless of inflation-state. Nevertheless, we note that the noise-like pattern observed in the parenchyma should arise due to the partially air-filled alveolar sacs. The  Fig. 1 Phase-contrast X-ray tomography of free-breathing mice using prospective respiratory gating. a Imaging arrangement, including anesthesia and gating equipment. Dashed lines with arrows indicate the direction of communication. The distances R 1 and R 2 are used to adjust the projection magnification and enables propagation-based phase-contrast X-ray imaging. b Zoom-in of the diaphragm (left) and ear (right) of a mouse captured in a 100 ms exposure projection (dark-and flat-field corrected). The phase-contrast produces edge-enhancement, seen here as a subtle bright contour around features (arrows). c Image acquisition flow with timings involved. Heartbeats can be observed as a low-amplitude modulation (~300-450 beats per minute (bpm)) of the measured respiratory movement. An adjustable delay (typically 300-500 ms) between respiration trigger and shutter opening ensures X-ray exposure (typically 100 or 150 ms) during the quiet respiratory phase. Graphs are raw data excerpts from a real experiment. d Histograms of pixel values in a representative dark-, flat-field and mouse projection acquired with the detector during a 100 ms exposure. Scale bars: 500 μm (b).
degradation of spatial resolution from uncompensated cardiac movement increases closer to the heart (cf. (A) in Fig. 4a), which limits the resolvability of the nearby smallest airways.
As we expected the tomographic spatial resolution to vary locally within the free-breathing mice, we performed a double CT acquisition (2 × 300 mGy) on a single animal by acquiring two projections per rotation angle (cf. Fig. 5a). By processing the two datasets separately, we acquired two unique tomographic reconstructions of the same mouse and timepoint. This allowed us to perform Fourier Shell Correlation (FSC) 34 on extracted regions-of-interests (ROIs) for estimating the local 3D resolution (cf. Fig. 5b). As each of the two reconstructed volumes contain residual motion blur averaged over the whole acquisition (in this case~20 min), real-domain comparisons of the volumes cannot directly yield any measure of motion blur and its influence on spatial resolution. However, as motion reduces the amount of high-frequency features in each of the volumes, FSC can shed light on the influence of residual motion on the local resolution by providing frequency-domain correlations. Specifically, we looked at three different regions (cf. example ROIs in Fig. 5c): (1) bones, where the internal structures should not be influenced by respiration (2) lung parenchyma including bronchioles, where we expect the most influence from respiration dynamics and (3) the 3D-printed support, which should be completely static during the scan. The FSC results are shown in Fig. 5d. In the support and bone regions, the 3D spatial resolution is predicted to be in the 30-μm-range (d min ), i.e., close to the system spatial resolution as measured earlier on the phantom. FSC of the lung parenchyma indicates that the intrapulmonary spatial resolution is somewhere between 34-100 μm (range where FSC coefficients approaches the threshold curve). This is in agreement with our direct measurements on large airways (cf. Fig. 4b) and confirms that we can reliably distinguish terminal bronchioles in the 60-100-μm-range in 3D. Contributing factors to the lower spatial resolution in the parenchyma should include cardiac and residual uncompensated respiratory motion.
Lastly, we extracted the tracheobronchial tree from the reconstructed volumes of the free-breathing mice using 3D segmentation methods. As seen in Fig. 6, the image quality is sufficient for segmenting sub-100-μm airways (animation of the 3D renderings shown in Supplementary Movie 3).

Discussion
Studies of animal lungs in their native state remain a key source of insight into human respiratory disease. We have shown that mechanical ventilation can be avoided for completely noninvasive and label-free imaging of the tracheobronchial tree in free-breathing mice using propagation-based phase-contrast X-ray tomography. Our laboratory arrangement provides imaging with 30-μm-range intrapulmonary spatial resolution which we have shown is sufficient for observing terminal bronchioles. Furthermore, the radiation dose delivered is below recommended levels for repeated exposure (500-800 mGy per week 25,26 ) and therefore allows for longitudinal imaging without adverse effects. By avoiding ventilator-induced lung-injury (VILI) this should find direct applications in studies of respiratory disease in mouse models, primarily those affecting the bronchioles (e.g., asthma, bronchitis, small airway disease, COVID-19). Even though we cannot resolve individual alveoli yet, studies of diseases presenting alveolar abnormalities as early manifestations (e.g., COPD, pneumonia) could benefit from the improved spatial resolution demonstrated here. Abnormalities on the alveolar level should result in macroscopic structural changes to the parenchyma which could be studied in greater detail with our proposed arrangement.
We have also shown that Fourier Shell Correlation (FSC) adapted from cryo-electron microscopy 34 and more recently super-resolution microscopy 35 is a powerful tool for estimating the local spatial resolution directly on in vivo tomographic data, which to our knowledge has been an unexplored tool in the in vivo imaging community. Here we further discuss mainly source-related limitations which need to be addressed for improving spatial resolution additionally towards alveolar-scale tomographic imaging of animal models in the laboratory.
Compared to previous studies, our results represent the highest in vivo tomographic resolution achieved in free-breathing mice  with a laboratory arrangement. Contrary to commercial preclinical micro-CT scanners which are often designed for fast lowdose scans and high-throughput studies, we designed our system for maximum spatial resolution and contrast in the lungs by using high-power liquid-metal-jet source technology enabling the long propagation distances required for phase-contrast imaging. We showed that our prospective gating successfully removes the bulk respiratory movement in free-breathing mice. Still, motion blur from heartbeats and residual diaphragm movement during the scan (cf. (B) in Fig. 4a) reduces our intrapulmonary spatial resolution by up to~3 × compared to the surrounding tissues and bones (cf. Fig. 5d). Prospective gating relies on the diaphragm and affected tissues relaxing to an identical position after each breath, which is an imperfect assumption when projections of freebreathing mice are acquired over hundreds of respiratory cycles.
Our results therefore indicate that we have reached the limit of what prospective respiratory gating alone can offer in terms of spatial resolution. Introducing cardiac gating would require 10ms-range exposure times (murine heart rates under anesthesia are in the range of~300-450 bpm) and has been important at synchrotrons for high-quality visualization of alveolar sacs (albeit at Gy-level radiation doses 12 ). In the laboratory setting, however, combining 10-ms-range exposure times with high resolution imaging is a challenge due to the limited brightness of microfocus X-ray sources (laboratory implementations of cardiac gating has not achieved sub-100-μm resolution 22,23,36 ). Smaller emission spots could additionally allow for higher geometric magnification resulting in increased spatial sampling in the object domain. An increase in photon flux will be necessary both for reducing exposure times and for preserving the per-pixel photon statistics at the increased spatial sampling, indicating that alveolar-level laboratory imaging using this approach is likely to come at the cost of higher dose (e.g., in the Gy-range, as observed at synchrotrons 12 ). Another way of reducing exposure times is to exploit binning in cases where the reconstructed voxels are fewfold smaller than the expected spatial resolution. Improved detectors based on photon-counting technology are likely to play an important role in keeping doses at acceptable levels by eliminating dark-counts (cf. Fig. 1d), read-out noise and improving quantum efficiency. The resolution of available photon-counting detectors should be high enough to resolve the edge-enhancement from the phase-contrast, but large pixel matrices are still expensive to manufacture as of today.
Another area of improvement lies with the phase-contrast strength, as increasing the latter would improve the air/tissue contrast and further facilitate segmentation tasks, as an example. In a laboratory setting, this typically translates into increasing both R 1 and R 2 and consequently the effective propagation distance z eff (¼ R 2 R 1 R 2 þR 1 ). However, increasing distances results in fewer X-ray photons reaching the detector at fixed exposure times. As photon statistics strongly influences image quality, one needs to either increase the exposure time per projection (not suitable for highquality gating) or increase the X-ray flux (preferred, but limited by available source technology). With the ongoing development of high-power microfocus X-ray sources like the recently released E1+ (Excillum AB, SE), we can envision future laboratory implementations with sources bright enough for increased z eff and thus stronger phase-contrast edge enhancement.
In summary, we have demonstrated that propagation-based phase-contrast X-ray tomography can be used to study the tracheobronchial tree down to the terminal bronchioles in freebreathing mice with radiation levels allowing for repeated imaging. As opposed to synchrotron implementations, this small-footprint arrangement (~1 × 2 × 3 m 3 ) should allow for broad adoption within the respiratory research community, specifically allowing for non-invasive longitudinal studies of murine lungs in their native state when forced breathing is not an option.

Methods
Animals. 5-10 week old female NMRI nude mice (BomTac:NMRI-Foxn nu , Taconic Biosciences, DK) weighing 20-30 g were used in this proof-of-principle study. Animals were housed in a temperature (22 ± 2 ∘ C) and humidity (55 ± 5%) controlled facility under a 12-h light/12-h dark cycle, with ad libitum access to food and water. Animal studies were carried out with approval from the regional animal ethics committee of Northern Stockholm, Sweden (protocol number 10579-2020), according to institutional and European guidelines for animal handling and research (EU Directive 2010/63/EU for animal experiments). The principles of the 3Rs were followed in our effort to minimize the number and suffering of animals.
Prior to imaging, mice were anesthetically induced in a chamber with 4-5% isoflurane mixed with either room-air or 100% O 2 from a concentrator (EverFlo, Philips, NL). No notable differences in anesthetic response and lung inflation were observed when comparing scans with either room-air or 100% O 2 as carrier. For anesthesia we used the Univentor 410 Anesthesia Unit and Univentor 2010 Scavenger Unit (Agnthos, SE). Following induction, mice were moved to a customdesigned 3D-printed vertical support for imaging and kept under maintenance anesthesia (0.8-1.5% isoflurane) using a facemask. The mice were wrapped in a thin layer of soft heat-isolating foam and secured to the support using non-sticking veterinary bandage, tight enough to minimize movement yet allowing adequate free-breathing. It is especially important to prevent image artifacts arising from vertical movement due to gravity during the course of the tomographic acquisition. We achieved this through our tight bandage wrapping, but other approaches could yield similar results. Body temperature was maintained during imaging using an on/off switchable infrared heating lamp. Opthalmic oinment (Oculentum simplex, APL, SE) was applied prior to imaging to prevent corneal drying. Post-imaging, animals were placed on a heating pad maintained at 40 ∘ C and monitored during anesthetic recovery. When respiration, mobility and reflexes returned to normal function the animals were returned to their cage.
Phase-contrast X-ray imaging. We used a liquid-metal jet X-ray source (D2, Excillum AB, SE) at 70 kV operating voltage and 200 W electron-beam power. The electron-beam spot was set to 20 × 80 μm resulting in an X-ray emission spot of 8.5 × 18 μm (FWHM). The X-ray beam was filtered using 1 mm Al to remove the strong emission peak at 9.2 keV (Ga K α ), resulting in a transmitted spectrum centered around the 24.2 keV peak (In K α ). X-ray exposure was controlled using a pneumatic shutter mounted before the filter. Source-to-object distance (R 1 ) was set to 27 cm and object-to-detector distance (R 2 ) to 135 cm, resulting in a geometric magnification of M ¼ R 1 þR 2 R 2 ¼ 6 onto the detector and an effective propagation The resulting X-ray flux on the detector plane was estimated to~3.5 × 10 5 X-ray photons/s/mm 2 (within the 15-30 keV spectral range). This corresponded to~90 X-ray photons (15-30 keV) per detector pixel during each 100 ms projection. Projection images were captured on a Xineos-2329 mammography CMOS flat panel detector (Teledyne DALSA, CA) with 49.5 μm pixels,~80 μm PSF and 23 × 29 cm 2 field of view (4608 × 5890 pixels). No pixel binning was used for the imaging experiments. The scintillator material of the detector was CsI with a non-specified thickness but most likely between 200-300 μm. In the object-plane, pixel sizes reduce to 8.25 μm due to the geometric magnification. During tomographic acquisition, mice are rotated 360 ∘ using a rotation stage (URS100 BCC, Newport, US) and projection images are acquired at 0. sampling used here is primarily limited by the suitable radiation dose per imaging session, but should otherwise be considered undersampled according to classical criteria (e.g., the Crowther criterion). To minimize both radiation dose and image distortion, as the Xineos-2329 lacks a global shutter, the detector acquisition was extended before and after the X-ray exposure controlled by the shutter, so that the detector was not exposed to X-rays during frame-cleaning/read-out. The shutter opening/closing time was measured using a 240 fps camera (iPhone 8, Apple Inc, US) and estimated to < 15 ms. Tomographic scan time was typically 12-20 min depending on the breathing-rate of the animal (i.e., due to respiratory-gating), but the actual X-ray exposure time was typically less than 2 min. An in-house developed graphical-user-interface (GUI) written in LabView (National Instruments, US) was used to control the tomographic acquisitions.
Prospective respiratory gating. We measured respiratory motion using a pneumatic pillow sensor (TSD110, BIOPAC Systems Inc., US) placed on the abdominal area of the animals and coupled to a monitoring system (MP160, BIOPAC Systems Inc., US) and gating unit (DTU200, BIOPAC Systems Inc., US) that triggers a TTL pulse at each breath. An in-house developed microprocessor-based (Arduino) controller uses the respiration trigger as input and outputs TTL pulses to the detector, shutter and rotation stage for synchronized respiratory-gated tomographic acquisition (cf. Fig. 1c). An adjustable delay following the trigger ensures X-ray exposure during the quiet respiratory phase. Gating settings (e.g., delay, exposure time) are transmitted to the controller using the same GUI used for tomographic acquisition.
Radiation dose estimates. We estimated the radiation dose delivered to mice during tomographic acquisition using the Monte Carlo simulation package MC-GPU 37 combined with the digital phantom DIGIMOUSE 38 . Using the X-ray source settings described above, the dose-rate was estimated to~4.2 mGy/s, calculated from the average whole-body dose (e.g., bones, muscle, adipose tissue) without tissue-specific weights (i.e., not converted to effective dose, measured in mSv). Deposited dose from Compton scattering and photoelectric absorption is taken into account. A tomographic acquisition with 720 projections, each with an X-ray exposure of 100 or 150 ms then results in an estimated total average dose of 300 or 450 mGy, respectively.
Image post-processing. Projection images were flatfield corrected to account for inhomogeneous X-ray illumination and interpolated to remove detector imperfections. For the respiratory motion assessment (cf. Fig. 2), the per-pixel mean square deviation (MSD) of a set of twenty consecutive projection images was calculated as MSDðm; nÞ ¼ ∑ p xðm; nÞ À x p ðm; where p denotes the projection index, m and n the pixel coordinates, x p the pixel intensity in projection p, and x the average pixel intensity in all projections. A 5 × 5 median filter was then applied to the calculated MSD maps with and without respiratory gating, and then scaled with the same factor for normalization. Here, MSD is a measure of pixel intensity deviation, i.e., a qualitative indicator of where motion is most present. In vivo local 3D resolution estimates. a Double computed tomography (CT) acquired on one animal (2 × 300 mGy), by acquiring two separate projection images for each rotation step (i.e., 2 × 720 = 1440 projections in total, each with 100 ms exposure). b The datasets were then processed individually (flat-field & dark-field correction, phase-retrieval and tomographic reconstruction). Fourier Shell Correlation (FSC) was then evaluated on regions-of-interests (ROIs) selected from both reconstructions. c Example of selected ROIs (2 × 2 × 2 mm 3 ) from bone, lung parenchyma and support regions. Grayscale in each ROI has been tuned separately for optimal contrast. d FSC as a function of the spatial frequency s (in line-pairs (lp) per mm) for each region plotted together with the corresponding half-bit threshold curve. Dashed lines represent raw FSC data for each ROI (here, number of ROIs (n) was 5 for each region). Solid lines with shaded regions represent the mean and standard deviation (std) of smoothed raw data. Scale bars: 250 μm (c).
Prior to tomographic reconstruction, corrected images were phase retrieved using Paganin's method 39 with the ratio δ/μ set to 4 × 10 −9 m for in vivo datasets and 5 × 10 −10 m for phantom datasets. The ratios were optimized for contrast and spatial resolution from the theoretical Δδ/Δμ values for the air/tissue interface in the lungs (~7 × 10 −9 m at 25 keV, i.e., source average energy) and air/Si interface in the phantom datasets (~1 × 10 −9 m). The two-material parameter is a valid extension of Paganin's method 40 . Tomographic reconstruction was performed with the FDK algorithm 41 using a commercial software (Octopus Reconstruction, Inside Matters NV, BE). Voxel sizes in the reconstructed volumes are 8.25 × 8.25 × 8.25 μm 3 , determined by the pixels in the projection images.
Spatial resolution analysis. We estimated the spatial-resolution of our imaging in both the real domain (e.g., line-profiles) and in the Fourier domain on both phantom data (Micro-CT Bar Pattern Air, QRM, DE) as well as directly on in vivo tomographic reconstructions. Fourier Ring Correlation (FRC) was performed on pairs of images using an ImageJ plugin 42 . Fourier Shell Correlation (FSC) was performed on pairs of ROIs (2 × 2 × 2 mm 3 or 242 × 242 × 242 voxels) using a Matlab (R2020a, Mathworks, US) implementation 43 modified for local resolution estimation (cf. Fig. 5). The latter was performed by multiplying the ROIs with a smooth-edged 3D Hann window function prior to FSC 44 . Representative ROIs (n = 5 for each region) were selected without spatial overlap. The raw FRC and FSC curves were smoothed using a 5-point moving average. We used the half-bit threshold as the criterion for estimating the spatial resolution for both FSC and FRC 34 .
3D segmentation and visualization. Prior to segmentation and 3D visualization, reconstructed volumes were binned and averaged by a factor of 2 into 16.5 × 16.5 × 16.5 μm 3 voxels and smoothed using a 3D median filter with 1 voxel radius for noise reduction while preserving feature edges. Ribs, lungs and tracheobronchial trees (cf. Fig. 6) were then segmented using a histogram-based region-growing segmentation method (Amira 6.3, Thermo Fisher Scientific, US). For airway segmentation, the region-growing parameters were tuned to the point where the smallest segmented features could still be distinguished from the surrounding lung parenchyma.
Reporting summary. Further information on research design is available in the Nature Research Reporting Summary linked to this article.

Data availability
A complete representative tomographic dataset (corresponding to Fig. 4) is available 45 . Any additional data will be made available upon request.