QUS Bug Report: Incorrect attenuation coefficient quantification & suggested fix on Clarius RF

[j1weeks]

---

name: General Issue and suggested Fix
about: Report bugs, request features, or ask questions
title: ''
labels: ''
assignees: ''

Issue Type

• Bug Report
• Feature Request
• Installation Issue
• Question/Help
• Documentation
• Other

Description

The current quantitative results for the attenuation coefficient on the Clarius system produce a very narrow range of values which are tightly clustered around the user-provided reference attenuation value. The source of the issue is an incorrect calculation of the axial depth, which flattens the slope $$p_0$$ in the following lines, located in qus/quantus/analysis/paramap/analysis_methods/attenuation_coef.py:

local_attenuation = ref_attenuation * f[f_idx] - (1 / 4) * p[0]
attenuation_coefficients.append(local_attenuation / f[f_idx])

The first line is equivalent to:

$$\mathrm{local_attenuation}(f)=\alpha_{\mathrm{ref}} f-\frac{1}{4} p_0$$

where

• $\alpha_{\text{ref}}$ is the reference attenuation (dB/cm/MHz)
• $f$ is the frequency (MHz)
• $p_0$ is the slope of the linear regression of spectral difference vs. depth

The reported attenuation estimate is:

$$\hat{\alpha}(f)=\alpha_{\text{ref}}-\frac{p_0}{4 f}$$

Final form of the estimator:

$$\hat{\alpha}(f)=\alpha_{\text{ref}}-\frac{p_0}{4 f}$$

In short, the actual estimator is just the user-provided reference attenuation value minus a correction term. This is broadly the same as the approach of Han et al. (https://doi.org/10.1002/jum.14537). Since f is fixed, the issue must be that $p_0$ is artificially small due to some issue earlier in the code.

Since $p_0$ is the slope of the log spectral difference versus depth, if the slope is small, the correction term will be even smaller. Therefore the estimate naturally collapses to the reference attenuation value.

Suggested Fix

In line 62 of attenuation_coef.py, the following line begins the scaling of the window depths:

# Convert window depths to cm
axial_res_cm = image_data.axial_res / 10

On the real Clarius C3 data, we find that image_data.axial_res = 1.994. But this is not a useful value for the remainder of the computation. If we examine the .yml file for the RF data, we find the following:

imaging depth: 99.70 mm
size: {samples per line: 1952,number of lines: 192,sample size: 2 bytes}

The axial spacing per RF sample, (dz), can be computed from the imaging depth and number of samples per line as:

$$ dz = \frac{\text{imaging depth (mm)}}{\text{samples per line}} $$

To convert to centimeters per sample:

$$ dz_\text{cm} = \frac{dz}{10} = \frac{\text{imaging depth (mm)}}{\text{samples per line} \times 10} $$

For example, if the imaging depth is 99.7 mm and there are 1952 samples per line as per our example:

$$ dz_\text{cm} = \frac{99.7}{1952 \times 10} \approx 0.00511 \text{ cm/sample} $$

Which is much smaller than the 0.1994 suggested by the current calculation. The current implementation is therefore inflating the RF_data_size by ~39x, resulting in an artifactually low slope.

The attenuation_coef function should use the .yml header data to determine the true resolution rather than using the image_data.axial_res, which is likely vestigial from some other metadata format. This will also be robust to changes in resolution due to different probes or presets.

I tested this fix locally; with hard-coded corrections to the resolution, a much more realistic range of values is achieved.

[j1weeks]

Addendum:

Please note that all information pertains to the latest main version and Clarius C3 data.

We discovered several discrepancies within the utc_loaders/clarius_rf/parser.py.

Swapped Axes in parser.py, reconstruct_bmode function

Above, we identified that axial_res was returning an incorrect value. This value is drawn from the scan-converted parser info:

        self.clarius_info_struct.yResRF =  self.clarius_info_struct.endDepth1*1000 / scBmodeStruct.scArr.shape[0]
        self.clarius_info_struct.axialRes = hCm1*10 / scBmodeStruct.scArr.shape[0]

However, print statements reveal that scArr.shape[0] = 500, while scArr.shape[1] = 859. Therefore the lateral and axial dimensions are swapped in the parser relative to our data.

Incorrect unit calculations in parser.py, reconstruct_bmode function

We investigated ~10 sample scans from different C3 probes over many days, all with ~10cm depth. Print statements revealed the following:

self.clarius_info_struct.endDepth1 * 1000: 997.0
scBModeStruct.scARr.shape[0]: 500
scBModeStruct.scARr.shape[1]: 859
hCm1 * 10 (supposedly depth): 820.7452057596531
self.clarius_info_struct.depth: 820.7452057596531
self.clarius_info_struct.axialRes: 1.6414904115193063

Therefore, the depth endDepth1 = 0.997, which means the depth is being measured in DECIMETERS. Multiplying this value by 1000 gives a depth of 99.7 meters, which is clearly nonsensical. We believe that the programmer believed the depth was being recorded in millimeters and multiplied by 1000 to give the depth in meters. The actual axialRes is calculated as

$$ self.clarius_/info_/struct.axialRes = hCm1 * 10 / scBmodeStruct.scArr.shape[0] $$

Where hCm1 (height in centimeters?) is also a nonsense value of 82.74 cm = 827.4 mm. The intent is further attested to by the following line and comment:

self.clarius_info_struct.depth = hCm1 * 10 #mm

which means that the depth calculation is also not working as intended and is very far off of the intended unit representation.

These issues are potentially carried through the rest of the code. However, the axialRes here is NOT the same as the image_data.axial_res mentioned in the issue above, which turns out to be

image_data.axial_res = 1.994 = 997/500 = self.clarius_info_struct.endDepth1 * 1000 / scBmodeStruct.scArr.shape[0]
= self.clarius_info_struct.yResRF

In summary, issues with image_data=axial_res are directly tied to the incorrect calculation of self.clarius_info_struct.yResRF in the parser, NOT the incorrect calculation of self.clarius_info_struct.axialRes. Depth data is recorded in decimeters, and subsequent unit conversions are done incorrectly. The shape axes also need to be transposed to correct for the swapping of the axes, but take note of the transposes already happening earlier in this function to ensure internal consistency.