MicroCT measurement

SH Simona Hankeova JS Jakub Salplachta TZ Tomas Zikmund MK Michaela Kavkova NH Noémi Van Hul AB Adam Brinek VS Veronika Smekalova JL Jakub Laznovsky FD Feven Dawit JJ Josef Jaros VB Vítězslav Bryja UL Urban Lendahl EE Ewa Ellis AN Antal Nemeth BF Björn Fischler EH Edouard Hannezo JK Jozef Kaiser EA Emma Rachel Andersson

This protocol is extracted from research article:

DUCT reveals architectural mechanisms contributing to bile duct recovery in a mouse model for Alagille syndrome

**
eLife**,
Feb 26, 2021;
DOI:
10.7554/eLife.60916

DUCT reveals architectural mechanisms contributing to bile duct recovery in a mouse model for Alagille syndrome

DOI:
10.7554/eLife.60916

Procedure

The system GE Phoenix v|tome|x L 240 (GE Sensing and Inspection Technologies GmbH, Germany) equipped with nanofocus X-ray tube (180 kV/15 W) was used for the tomographic measurements that were carried out in the air-conditioned cabinet (fixed temperature 21°C). The samples were adapted for this temperature before the measurement to prevent any thermal expansion effect. To prevent any sample motion during the scanning, the samples were placed in 15 ml Falcon tube, filled with 1% agarose gel. The tomographic reconstruction of acquired data was performed using GE phoenix datos|x 2.0 software. The voxel resolution was fixed for all the adult liver samples at 12 µm, except one (sample #2401, 8 µm). For all the P15 liver samples the voxel resolution was fixed at 6.5 µm and for the lung lobe sample at 8 µm. Detailed overview of used acquisition parameters is stated in Table 2.

*Flat panel DXR250 (2048 px ×2048 px, pixel size 200 μm).

^{†}Flat panel dynamic 41|100 (4048 px ×4048 px, pixel size 100 μm with binning 2).

The identification and segmentation of both tubular systems (e.g. bile duct (BD) and portal vein (PV) for liver samples and pulmonary artery and airways for lung samples) in each CT cross-section was necessary for further analysis of each system. The segmentation was based on differential contrast between the resin and the soft tissue. Two different resins were used to identify the individual tubular systems. The differential contrast for system identification was highly dependent on the freshness of the MICROFIL.

The resin segmentation was performed by global thresholding in VG Studio MAX 3.3 (Volume Graphics GmbH, Germany) software together with manual corrections where necessary. The threshold value was determined based on the histogram shape and visual evaluation of a selected cross-section (Figure 8A,A'). The resin cast, especially when using an old MICROFIL, could contain artefacts caused by poor contrast, insufficient filling or high injection pressure. These artefacts include air bubbles in resin, non-homogenous resin contrast (caused by mixing of blue and yellow MICROFIL) and resin leakage due to lumen rupture (probably caused by high pressure during injection), (Figure 8B). Therefore, the thresholding step was supplemented with manual corrections to create smooth, continuous and solid canal masks. Furthermore, the cut-off for the smallest distal canal included in the mask is considered an area of at least four voxels.

Next, the individual tubular systems were identified in the segmented resin mask (Figure 1—figure supplement 12D and G). When fresh MICROFIL was used, it was possible to identify each system by global thresholding with threshold value determined based on histogram shape and visual evaluation of CT data. However, in case of the old MICROFIL (used for adult liver samples) the resin absorption properties were not distinguishable in CT data. In most regions, both tubular systems blended with each other in one continuous region (Figure 8C, D left panel). Manual segmentation was therefore necessary to ensure the correct identification of both systems (Figure 8C, D right panel). The manual segmentation was performed by outlining the BD regions in every slice of the CT data. VG Studio automatically creates 3D render based on the regions outlined in CT sections.

The adult liver was analyzed using a custom-written algorithm and freely available Matlab codes (Version R2017a, The MathWorks Inc, Natick, MA). The algorithm was designed to analyze morphological parameters of the BD and PV systems, and is compatible with the 3D binary masks. Two separate masks of BD and PV system were generated and the analysis was divided in two independent parts. First, the analysis of the entire portal vein and biliary system and second, analysis of the corresponding main branch (=the longest branch) of each system (Figure 6—figure supplement 1A). For detailed analysis and comparison of the whole system versus only the main branch, two algorithms, described by the diagrams in Figure 8F, G, were developed. They differ in the input data and the evaluated parameters. For both algorithms, the first step is to create a 3D skeleton of the input binary mask.

The 3D skeleton was derived using the homotopic thinning algorithm described in Lee et al., 1994 specifically optimized for Matlab implementation by Kollmannsberger (Kerschnitzki et al., 2013), (online source: https://www.mathworks.com/matlabcentral/fileexchange/43400-skeleton3d). Calculated 3D medial axis skeleton was subsequently converted to a network graph, using algorithm described in Kerschnitzki et al., 2013 (online source: https://www.mathworks.com/matlabcentral/fileexchange/43527-skel2graph-3d). Resulting network graph is formed by nodes and links between them (Figure 1—figure supplement 12D and G).

Liver was separated into three regions: R1 (approximately hilum), R2 (approx. intermediate) and R3 (approx. periphery). The optimal region size was calculated as a summary of average PV branch lengths per generation and divided by three (detailed branching generation subdivision Table 3 for PV and Table 4 for BD). Subsequently the branching generations were assigned to a region by matching the summary of branch length per generation to optimal region size. The same region was applied for both PV and BD analysis with the exception of one *Jag1 ^{Ndr/Ndr}* sample (#1 or 2714) where the bile duct optimal region size was greater than portal vein.

Regions were assigned in order that total average lengths of each generation within each region yielded an equal size of R1, R2, and R3. Each sub-column represents one animal.

The distribution of BD branches within each region were quantified based on regions defined by the PV system (Table 4). Each sub-column represents one animal.

Branching points analysis was programmed in Matlab to analyze the distance between BD and PV branching point. This parameter was calculated using 3D Euclidean distances between the BD branching points and the nearest branching point from PV system. The data is represented as cumulative sum of percentage of BD branching point at a given distances between BD and PV branching points (from 0.015 mm to 3 mm).

For branch length analysis the structure of BD and PV trees were first reconstructed in 3D, using the Analyze Skeleton toolbox in ImageJ, which provided the three-dimensional coordinates of all branch points for both BD and PV, as well as the connectivity of the graph. Next, we computed for each branch the length along its path to the Euclidean distance between its extremities (branch points). To calculate the generation number of branches (both for the PV and BD structures), we manually defined the origin of the ducts and vessels as generation 1, and computed generation number as the number of generation branches separating a given branch from the origin. To distinguish side branching events, we calculated the angle between a branch and its ‘parent’ by computing the dot product p of both their unit vectors. A branch with p>0.95 with its parent branch was considered to belong to the same generation. We then computed distributions of length for the BD and PV structures as a function of generation number. Each generation was assigned to a region R1, R2, or R3.

The branch distribution in each region was calculated as a summary of number of branches per generation in a given region. The sum of branch numbers of each region is displayed as a proportion of the total number of branches per sample.

The number of bi-furcations (i.e. one input and two outputs), tri-furcations (i.e. one input and three outputs) and quadri- and more-furcations (one input and more than three outputs), were assessed in Matlab based on binary mask skeleton nodes that were divided into endpoints and branching points. Branching points closer than 0.2 mm (this threshold value was derived based on visual assessment and knowledge of the system) were merged together and further represented by one node.

To evaluate the gap between BD and PV the surface distances were calculated in Matlab for each BD skeleton point by detecting the nearest PV skeleton point and connecting the two points with a line and measuring the non-resin area on this line (zero area in the input binary masks). Surface distance was then calculated using 3D Euclidean distance between the detected non-resin voxel coordinates. The data is represented as cumulative sum of percentage of BD at a given distance from PV (from 0.015 mm to 1.5 mm). The maximum distance between BD and PV for each liver sample was depicted in a separate graph.

To quantify length and tortuosity, total (curved) and theoretical (chord) lengths were measured in Matlab for the whole system length and for the corresponding main branch. The curved length was defined as a cumulative sum of 3D Euclidean distances between neighboring graph points (i.e. links forming points) multiplied by voxel size. The chord length was defined as cumulative sum of 3D Euclidean distances between neighboring nodes multiplied by voxel size. The chord length therefore reflects system length where any nodes are connected by links with the shortest possible length. To analyze the relationship between BD and PV a length of the BD was divided by a PV length (curved or chord). Tortuosity was calculated as curved length divided by chord length of the same system and distributed into regions based on the generation number as previously described. Tortuosity was assessed in %, as BD and PV are not straight lines the actual tortuosity measurements were subtracted by 100% (perfectly straight line).

Total system volume was calculated in Matlab by multiplying a number of voxels representing PV or BD by volume of one voxel. The relationship between BD and PV volumes was addressed by dividing the BD volume by PV volume.

The main branch diameter was calculated in Matlab every 1.5 mm along the total length of the main branch. The radius was defined as the minimal distance from the skeleton to the segmented area boundary in the input binary mask (i.e. border between background and area of interest). This boundary was calculated using a two-step procedure. In the first step, the input map was eroded using a 3D spherical shaped structural element with one pixel radius. Subsequently, the eroded area was subtracted from the original binary mask. This resulted in a binary mask representing the boundary between the background and the area of interest. One radius value at a given skeleton point was then expressed as the minimum distance from that point to the mask boundary. This was calculated, using the minimal value search in the intersection of the boundary mask and the distance map from that point. The distance map from a given skeleton point was calculated as 3D Euclidean distance of the spatial coordinates. Subsequently the diameter value was calculated as the minimum distance to the boundary area multiplied by 2. To avoid any misrepresentation, the one final diameter value at a given point (every 1.5 mm of branch length) was calculated as a mean value of a diameter at that point and diameters at four neighboring points (two on each side). PV and BD diameters were divided into three areas: hilum, intermediate and periphery. Hilar region represents distance from 0 to 1.5 (sample #2401) or 0–3 mm (other samples), Intermediate region is from 3 to 6 mm (sample #2401) or 4.5 mm – 9 mm (other samples), Periphery is from 7.5 mm – 9 mm (sample #2401) or 10.5–13.5 mm (other samples). BD to PV diameter ratio was calculated by dividing BD diameter at a given region by PV dimeter of the same region.

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.

Note: The content above has been extracted from a research article, so it may not display correctly.

Q&A

Your question will be posted on the Bio-101 website. We will send your questions to the authors of this protocol and Bio-protocol community members who are experienced with this method. you will be informed using the email address associated with your Bio-protocol account.