Quantitative Assessment of the Echogenicity of a Breast Tumor Predicts the Response to Neoadjuvant Chemotherapy

Simple Summary B-mode US is a widely available, inexpensive, and non-invasive technique. This method is used in monitoring the neoadjuvant chemotherapy (NAC) in breast cancer (BC). In the presented study we combined the result from B-mode ultrasound examination with quantitative information about the characteristics and structure of the tissue in predicting the response to neoadjuvant chemotherapy in BC patients. We used echogenicity (ΔEcho) as B-mode features and the Kullback-Leibler divergence (ΔKLD) method as a quantitative parameter to provide information on changes in image echogenicity, to determine differences between the distributions of the ultrasound echo amplitude from tumor during NAC. The ΔKLD parameter alone is an accurate predictor of response to treatment after the second course of therapy (cut-off ≥70%, AUC = 0.85). Combining both parameters (ΔKLD and ΔEcho) led to an increase in sensitivity without significant deterioration of other statistical parameters and allowed to accurately predict non-responding tumors. Abstract The aim of the study was to improve monitoring the treatment response in breast cancer patients undergoing neoadjuvant chemotherapy (NAC). The IRB approved this prospective study. Ultrasound examinations were performed prior to treatment and 7 days after four consecutive NAC cycles. Residual malignant cell (RMC) measurement at surgery was the standard of reference. Alteration in B-mode ultrasound (tumor echogenicity and volume) and the Kullback-Leibler divergence (kld), as a quantitative measure of amplitude difference, were used. Correlations of these parameters with RMC were assessed and Receiver Operating Characteristic curve (ROC) analysis was performed. Thirty-nine patients (mean age 57 y.) with 50 tumors were included. There was a significant correlation between RMC and changes in quantitative parameters (KLD) after the second, third and fourth course of NAC, and alteration in echogenicity after the third and fourth course. Multivariate analysis of the echogenicity and KLD after the third NAC course revealed a sensitivity of 91%, specificity of 92%, PPV = 77%, NPV = 97%, accuracy = 91%, and AUC of 0.92 for non-responding tumors (RMC ≥ 70%). In conclusion, monitoring the echogenicity and KLD parameters made it possible to accurately predict the treatment response from the second course of NAC.


Introduction
Breast cancer (BC) is a disease of significant social importance and is the most common malignant neoplasm in women in Poland and worldwide. The incidence of BC in women Experiments using ultrasound data collected from cell pellets exposed to chemotherapeutics that induce apoptosis have shown a significant effect of apoptosis on ultrasound scattering and the echogenicity of ultrasound images [10,11]. These results encourage the use of Quantitative Ultrasound (QUS) in the context of monitoring the effects of chemotherapy.
The possibilities of using QUS methods in NAC monitoring have recently been widely explored. For example, Quiaoit et al. proposed a multi-parameter approach to the problem of detecting a positive response to therapy [12]. In this study, the authors generated parametric images using the spectral and backscattering parameters of the signals scattered within the tumor. Then they used the texture features of these images to develop multiparametric models to predict pathological tumor response.
The best model was able to predict positive response to NAC with area under the ROC curve (AUC) value equal 0.87 after the first week of the therapy. Dasgupta et al. used textural features of QUS parametric images of spectral and scattering parameters of a tumor for predicting the response to NAC with AUC of 0.86 [13]. In other work, Tadayyon et al. showed the potential of combining QUS analysis with artificial neural networks as a source of diagnostic biomarkers [14]. The authors reported high performance of their model in predicting the pathological response (AUC = 0.96) and 5-year recurrence-free survival (RFS) of patients (AUC = 0.89) prior to the start of treatment. Despite high efficiency reported in those studies, there is still much room for improvement. Each new parameter that is reported to perform decent in predicting the response to NAC, is valuable, as it can potentially be included in a multi-parametric model, improving its performance.
In our study, we focused on the characterization of the tumor echogenicity, which is related to its structure and cellularity, as is the ADC in the MRI image [15,16]. In the classic B-mode examination, echogenicity is assessed in relation to adipose tissue, which makes the assessment more objective and independent of the ultrasonic scanner settings. The echogenicity of a given area of the examined tissue translates into the amplitude of the scattered signal and affects the average value of the amplitude distribution. Tissue signal amplitude distribution depends on tissue structure, and their mean values may differ depending on the tissue. Moreover, with a similar mean value, they can describe different tissue structures. Therefore, the differences between the distributions are more accurately determined using statistical measures of variation than by their mean values.
We used the Kullback-Leibler divergence (kld) method [17] to determine differences between the distributions of the ultrasound echo amplitude from tumor tissue undergoing successive cycles of NAC therapy, and, ultimately, to predict the effects of NAC. Amplitude distributions were determined from raw radio-frequency (RF) signals, which ensured independence from the processing methods used by scanners for image enhancement. At the same time, we investigated the possibility of predicting the effects of NAC using changes in B-mode image echogenicity and tumor size and the combination of both parameters. The results were assessed in relation to postoperative residual malignant cells (RMC) rate.

Patients
This prospective, single-center study was conducted in accordance with the Declaration of Helsinki, and the protocol was approved by the Ethics Committee of Maria Skłodowska-Curie National Institute of Oncology, Scientific Centre, Warsaw, Poland (Project identification code 49/2018). All participants provided informed consent for inclusion before participation in the study.
Breast ultrasound examinations were performed on 39 patients with a total of 50 BCs (seven women had bifocal lesions and two had trifocal lesions) from April 2016 to July 2020.
The inclusion criteria were as follows: inclusion in the NAC by an oncologist during multidisciplinary meeting (MDM); maximum tumor diameter <4 cm, minimum tumor diameter 5 mm; and multicenter ≤3 if in another quadrant and/or breast, regardless of the immunohistological subtype and lymph node status. Ultrasound examinations were performed before NAC and 7 days after subsequent NAC courses (first to fourth). NAC was administered according to international guidelines, according to a previously detailed protocol [18]. AC (doxorubicin, cyclophosphamide) was used from the first to the fourth course, and then treatment was continued with taxol. Patients with HER2+ receptor-positive tumors were treated with trastuzumab and taxol. One patient with a history of contralateral BC 5 years prior to the study, was treated with AT (doxorubicin and docetaxel), 37 patients underwent mastectomy, 2 underwent BCT, 37 underwent lymphadenectomy, and 2 underwent sentinel lymph node dissection (SLND) at the end of chemotherapy. Four patients did not complete the study, which resulted in modification of the NAC data. All remaining patients underwent a full course of NAC therapy.

Histology
In all patients, core needle biopsies (CNB) were performed after administration of 2% lidocaine before qualification to NAC treatment (14GA diameter biopsy needle, three to five cores; Pro-Mag). An experienced pathologist with 25 years' experience in breast oncological pathology assessed the cores from CNB and excisable tumors or residual intramammary target lesions with clip markers after NAC. Cancer subtypes (molecular subtypes and grade of malignancy) were obtained by pathological assessment after CNB. RMC rates in the residual malignant tumor were estimated by a pathologist based on postoperative histopathology and served as a reference of tumor response to NAC. RMC is one of the six features of the residual cancer burden (RCB) index that is a quantitative parameter evaluating the response to NAC that can be calculated using an on-line calculator. The percentage of RMC ranged from 0 to 100% and was determined after surgery [19]. In the statistical analysis, we evaluated the following cut-off values for RMC: ≤30% for responding tumors and ≥70% for non-responding tumors.

Ultrasonic Data Registration
The acquisition of RF ultrasound echoes from patients was performed using an Ultrasonix SonixTOUCH ultrasound scanner (Ultrasonix Medical Corporation, Richmond, BC, Canada). The breast lesions were assessed according to American College of Radiology (BI-RADS-lexicon) and the standards of Polish Ultrasound Society [20,21]. The scanner was equipped with an ultrasound research interface enabling the recording of raw postbeamformed RF data. The measurements were made using a linear probe (L14-5/38) at a frequency of 10 MHz. In accordance with the measurement protocol, the ultrasound RF data were recorded for each tumor in the following planes: radial, radial +45 • , anti-radial, and anti-radial +45 • . Classic amplitude images (B-mode) were generated from these data, with no additional post-processing applied. This is important because commercial scanners usually perform additional image processing to improve the images, for example, to increase contrast. Such operations influence the statistics of the image amplitudes to such a degree that they are difficult to assess. Therefore, in these studies, non-post-processed data were used.
Data registration was performed before the start of treatment and 7 days after each course of chemotherapy. In each B-mode image, a tumor was indicated by a radiologist with 20 years of experience in breast ultrasound. Data from four cross sections of the tumor were analyzed as a single set. This approach allowed for a more detailed analysis of changes in the tumor microstructure than was possible by examining only a single cross-section.

Quantitative Analysis of Ultrasound Data
The physical properties and spatial distribution of scatterers influence the statistical distribution of the amplitudes of the received echoes. A diagram of the ultrasonic signal generation is shown in Figure 1. Thus, the amplitude distribution of the ultrasonic signal received by the transducer contains information about the microstructure of the examined tissue. For example, if breast microcalcifications are present (a source of strong ultrasound scattering), higher amplitude values may be expected (Figure 1). changes in the tumor microstructure than was possible by examining only a single crosssection.

Quantitative Analysis of Ultrasound Data
The physical properties and spatial distribution of scatterers influence the statistical distribution of the amplitudes of the received echoes. A diagram of the ultrasonic signal generation is shown in Figure 1. Thus, the amplitude distribution of the ultrasonic signal received by the transducer contains information about the microstructure of the examined tissue. For example, if breast microcalcifications are present (a source of strong ultrasound scattering), higher amplitude values may be expected (Figure 1). Changes in the microstructure of the tumor occur as a result of the NAC effect. Remodeling is a multi-stage process, involving changes to cancer cells (cellular organization, cell death) and changes in the stroma. The degree of change in tumor depends on individual factors. However, large changes in the microstructure of responding tumors and small changes in tumors resistant to treatment with NAC can be expected. Therefore, it can be assumed that the degree of change in the amplitude distributions of the received signals is directly related to the degree of change at the microstructure level. This assumption forms the basis of the method presented in this study. The amplitude distributions were estimated using kernel density estimation (KDE) [22,23]. Amplitude distributions determined for responder (RMC ≤ 30%) and non-responder (RMC ≥ 70%) patients are shown in Figure 2, respectively. Changes in the microstructure of the tumor occur as a result of the NAC effect. Remodeling is a multi-stage process, involving changes to cancer cells (cellular organization, cell death) and changes in the stroma. The degree of change in tumor depends on individual factors. However, large changes in the microstructure of responding tumors and small changes in tumors resistant to treatment with NAC can be expected. Therefore, it can be assumed that the degree of change in the amplitude distributions of the received signals is directly related to the degree of change at the microstructure level. This assumption forms the basis of the method presented in this study. The amplitude distributions were estimated using kernel density estimation (KDE) [22,23]. Amplitude distributions determined for responder (RMC ≤ 30%) and non-responder (RMC ≥ 70%) patients are shown in Figure 2, respectively. Average probability density (PDF) function for backscatter amplitude obtained for responders (RMC ≤ 30%) and non-responders (RMC ≥ 70%). Data were recorded before chemotherapy (NAC0) and after cycle 3 (NAC3).
The Kullback-Leibler divergence ( ) [13] was used to quantify the differences between two distributions of amplitudes (where hn and hm denote the amplitude distributions after the n-th and m-th treatment courses), which is given by the equation:  The Kullback-Leibler divergence (kld) [13] was used to quantify the differences between two distributions of amplitudes (where h n and h m denote the amplitude distributions after the n-th and m-th treatment courses), which is given by the equation: where i is a sample index in the amplitude distribution. The values of kld are always non-negative; a value of 0 indicates that the distributions h n and h m are identical. The kld was used to define three parameters: where h is the estimated amplitude distribution, its subscript of 0 indicates pre-treatment data, and a positive subscript n indicates the NAC course at which the data were recorded. The use of logarithms in the above formulas results from highly asymmetric distributions of the KLD parameters if they are devoid of logarithm. After application of the logarithm, these distributions were close to normal.
To evaluate the effectiveness of NAC therapy, three quantitative parameters, all based on the KLD method, were used. The KLD 0 parameter estimates the changes in the amplitude distributions in the tumor after particular NAC courses in relation to the distributions before the start of therapy. The KLD 1 parameter describes the changes in amplitude distributions with respect to the amplitude distribution after the first NAC, and the ∆KLD parameter describes the differences between the distributions in relation to the difference between the distribution after the first NAC and the distribution before treatment initiation.

Tumor Echogenicity
Tumor echogenicity (Echo) was assessed by a physician based on gray-level standard Bmode images compared to fat tissue in the preglandular zones. The following echogenicity levels were assigned to the images of each tumor: 1. Hypoechoic 2.

Isoechoic 4. Hyperechoic
The Echo change (∆Echo) was determined in relation to pre-treatment value Echo(0) using the following equation: where n indicates the NAC course

Tumor Volume
The volume (V) of the tumor was calculated assuming an ellipsoidal shape using the following equation: where w, d, and h represent the width, depth, and height of the tumor, respectively, assessed based on the B-mode images. The volume change (∆V) was determined in relation to pretreatment value V(0) using the following formula: where n indicates the NAC course. For the same reason as KLD, ∆V is subjected to logarithm. In order to enable an intuitive evaluation of volume changes, a percentage change in volume ∆V % was defined as:

Statistical Analysis
All of the considered parameters were examined in terms of their correlation with the RMC at particular stages of therapy. As none of the abovementioned parameters met the normality criterion at all stages of the therapy (Shapiro-Wilk test, significance level 0.05), the Spearman's rank correlation coefficient r S was used to assess the correlation. Subsequently, individual parameters were analyzed in the context of tumor classification as "responding" and "non-responding." Tumors with RMC ≤ 30% were considered responders, and those with RMC ≥ 70% were considered non-responders. Each of the above parameters was analyzed as a single-parameter classifier. In addition, the ∆Echo and ∆V parameters were combined in a two-parameter classifier using linear discriminant analysis (LDA), which is a generalization of Fisher's linear discriminant [24]. The tested classifiers were cross-validated using the leave-one-out (LOO) method [25]. The only exception was the classifier based solely on ∆Echo, where, because of the coarse discretization of the parameter, the cross-validation algorithm led to a distortion of the analysis results and was therefore omitted.
The effectiveness of the classification was assessed based on the receiver operating characteristics (ROC) curve, which demonstrated a trade-off between the true positive rate (TPR) and the false positive rate (FPR) for a classifier. The area under the ROC curve (AUC) was used for an overall assessment of the entire ROC curve [26]. For a random classifier, the AUC is close to 0.5, whereas for the ideal classifier, the AUC is equal to 1. The confidence intervals for the AUC were determined using the bootstrap method. The number of bootstrap samples was 10 4 , and the confidence level was 0.95. Each classifier was assessed in terms of its sensitivity, specificity, accuracy, positive predictive value (PPV), and negative predictive value (NPV). These features were determined for the optimal operating point, which was chosen as the point on the classifier's ROC curve closest to the point (FPR = 0, TPR = 1) in the Euclidean sense [27].

Tumor Size
The mean size of the lesion before the treatment (n = 50) was 5.0 cm 3 (median, 2.8 cm 3 ; range, 0.04-27 cm 3 ; SD, 5.7 cm 3 ), and after three courses of NAC (n = 47) it was 2.0 cm 3 (median, 0.67 cm 3 ; range, 0.02-14 cm 3 ; SD, 2.8 cm 3 ). When analyzing changes in the volume of neoplastic tumors, decreases were observed in the entire group of tumors after the first and second courses of NAC, while a reversal of this trend was observed after the third course of NAC. For tumors that responded poorly to treatment (RMC ≥ 70%), an increase in tumor volume was observed in relation to the previous courses, although they were still smaller than the pre-treatment measurement ( Figure 5).
(median, 0.67 cm 3 ; range, 0.02-14 cm 3 ; SD, 2.8 cm 3 ). When analyzing changes in the volume of neoplastic tumors, decreases were observed in the entire group of tumors after the first and second courses of NAC, while a reversal of this trend was observed after the third course of NAC. For tumors that responded poorly to treatment (RMC ≥ 70%), an increase in tumor volume was observed in relation to the previous courses, although they were still smaller than the pre-treatment measurement ( Figure 5).

The Kullback-Leibler Divergence Based Parameters
The statistics of the Kullback-Leibler divergence (kld)-related parameters are shown in Figure 6. Distributions of each parameter show improvement in the separation of responders and non-responders with subsequent courses of NAC. The best separation is observed for the ΔKLD parameter.

The Kullback-Leibler Divergence Based Parameters
The statistics of the Kullback-Leibler divergence (kld)-related parameters are shown in Figure 6. Distributions of each parameter show improvement in the separation of responders and non-responders with subsequent courses of NAC. The best separation is observed for the ∆KLD parameter.

Correlation
Prior to analyzing the usefulness of the parameters as predictors of the efficacy of NAC, their correlations with the RMC results were assessed. After the first course of NAC, the ΔEcho, ΔV, and KLD0 parameters showed weak and mostly statistically insignificant correlations, which only slightly improved with subsequent NAC courses. The KLD1 and ΔKLD parameters, which were available starting from the second course of NAC, demonstrated statistically significant, mostly moderate or strong correlations (Figure 7).

Correlation
Prior to analyzing the usefulness of the parameters as predictors of the efficacy of NAC, their correlations with the RMC results were assessed. After the first course of NAC, the ∆Echo, ∆V, and KLD 0 parameters showed weak and mostly statistically insignificant correlations, which only slightly improved with subsequent NAC courses. The KLD 1 and ∆KLD parameters, which were available starting from the second course of NAC, demonstrated statistically significant, mostly moderate or strong correlations (Figure 7).

Correlation
Prior to analyzing the usefulness of the parameters as predictors of the efficacy of NAC, their correlations with the RMC results were assessed. After the first course of NAC, the ΔEcho, ΔV, and KLD0 parameters showed weak and mostly statistically insignificant correlations, which only slightly improved with subsequent NAC courses. The KLD1 and ΔKLD parameters, which were available starting from the second course of NAC, demonstrated statistically significant, mostly moderate or strong correlations (Figure 7).  A table demonstrating the resulting Spearman correlations and p-values is available in Appendix A (Table A1).

Classification
Two RMC thresholds (RMC ≤ 30% and RMC ≥ 70%) were analyzed in order to divide the tumors into "responder" and "non-responder" groups. Figure 8 shows the AUC values together with confidence intervals for the analyzed parameters, as a function of chemotherapy courses. The corresponding ROC curves are included in the Appendix A ( Figure A1). The Appendix A also includes complete tables of classification performance of responders and non-responders (Tables A2 and A3, respectively). Their essential parts are presented in the main body of the article as Tables 3-5, for later discussion. Table 3. Results of statistical analysis of changes in echogenicity, volume, and ∆KLD parameter after two, three, and four courses of NAC for non-responding tumors (RMC ≥ 70%).   Of the B-mode parameters, assessed separately and in combination, the assessment of changes in echogenicity (ΔEcho) most accurately predicted the group of non-responding tumors with RMC ≥ 70% (Table 3) after the third NAC course (AUC = 0.81), with a sensitivity of 64% and specificity of 86%. To evaluate the volume change (ΔV), lower AUC values were obtained after both the third (AUC = 0.63) and fourth (AUC = 0.79) courses for RMC ≥ 70%.

Measure of Performance
We found that there was a slight improvement in accuracy, specificity, and PPV with the combined model of both parameters (change in echogenicity and volume) (accuracy 87% vs. 80%, specificity 94% vs. 86%, 58% vs. 78%), with no improvement in sensitivity (64%) ( Table 5). The results for all of the RMC cutoff values are shown in Appendix A (Tables A2 and A3).
Starting from the second course of chemotherapy, the quantitative parameter ΔKLD made it possible to predict the response to treatment with AUC ≥ 0.85. The ROC curves obtained after the second, third, and fourth courses of chemotherapy are presented in Figure 9. As in the case of echogenicity, the ΔKLD brought the best results after the third treatment course (AUC = 0.90; sensitivity of 82%; specificity of 94%). All of the classifiers characteristics are listed in Table 3.  Of the B-mode parameters, assessed separately and in combination, the assessment of changes in echogenicity (∆Echo) most accurately predicted the group of non-responding tumors with RMC ≥ 70% (Table 3) after the third NAC course (AUC = 0.81), with a sensitivity of 64% and specificity of 86%. To evaluate the volume change (∆V), lower AUC values were obtained after both the third (AUC = 0.63) and fourth (AUC = 0.79) courses for RMC ≥ 70%.
We found that there was a slight improvement in accuracy, specificity, and PPV with the combined model of both parameters (change in echogenicity and volume) (accuracy 87% vs. 80%, specificity 94% vs. 86%, 58% vs. 78%), with no improvement in sensitivity (64%) ( Table 5). The results for all of the RMC cutoff values are shown in Appendix A (Tables A2 and A3).
Starting from the second course of chemotherapy, the quantitative parameter ∆KLD made it possible to predict the response to treatment with AUC ≥ 0.85. The ROC curves obtained after the second, third, and fourth courses of chemotherapy are presented in Figure 9. As in the case of echogenicity, the ∆KLD brought the best results after the third treatment course (AUC = 0.90; sensitivity of 82%; specificity of 94%). All of the classifiers characteristics are listed in Table 3. The influence of the selection of the RMC threshold on the efficiency of the ΔKLDbased classifier is shown in Figure 10. The AUC increases with the adopted RMC threshold, which means that non-responding tumors can be indicated more efficiently. Detection of the tumors with the worst response to NAC (RMC ≥ 90%) is characterized with AUC significantly exceeding 0.90. While the classification of responders with use of the ΔKLD parameter appears more difficult, it still allows to achieve AUC = 0.84 for RMC ≤ 30% (Table 4). Figure 10. AUC of the ΔKLD-based classifier as a function of the adopted RMC threshold level.

Discussion
In the current study, ultrasound images were analyzed to monitor the response to NAC after the first four courses. Two types of images were included: traditional B-mode images and amplitude images generated from reconstructed (post-beamformed) RF data. The influence of the selection of the RMC threshold on the efficiency of the ∆KLDbased classifier is shown in Figure 10. The AUC increases with the adopted RMC threshold, which means that non-responding tumors can be indicated more efficiently. Detection of the tumors with the worst response to NAC (RMC ≥ 90%) is characterized with AUC significantly exceeding 0.90. While the classification of responders with use of the ∆KLD parameter appears more difficult, it still allows to achieve AUC = 0.84 for RMC ≤ 30% (Table 4). The influence of the selection of the RMC threshold on the efficiency of the ΔKLDbased classifier is shown in Figure 10. The AUC increases with the adopted RMC threshold, which means that non-responding tumors can be indicated more efficiently. Detection of the tumors with the worst response to NAC (RMC ≥ 90%) is characterized with AUC significantly exceeding 0.90. While the classification of responders with use of the ΔKLD parameter appears more difficult, it still allows to achieve AUC = 0.84 for RMC ≤ 30% (Table 4). Figure 10. AUC of the ΔKLD-based classifier as a function of the adopted RMC threshold level.

Discussion
In the current study, ultrasound images were analyzed to monitor the response to NAC after the first four courses. Two types of images were included: traditional B-mode images and amplitude images generated from reconstructed (post-beamformed) RF data. The purpose of using the RF data was to avoid the unknown impact of image post-pro- Figure 10. AUC of the ∆KLD-based classifier as a function of the adopted RMC threshold level.

Discussion
In the current study, ultrasound images were analyzed to monitor the response to NAC after the first four courses. Two types of images were included: traditional B-mode images and amplitude images generated from reconstructed (post-beamformed) RF data. The purpose of using the RF data was to avoid the unknown impact of image post-processing algorithms implemented in commercial scanners.
In traditional B-mode images, changes in volume and echogenicity were analyzed relative to the baseline study before treatment. Quantitative ultrasonic examination was based on the KLD analysis. Three parameters, KLD 0 , KLD 1 , and ∆KLD, were used to assess the discrepancy between the distributions of echogenicity of the images after particular NAC courses compared to the distributions for amplitude images determined: (1) before the start of treatment and (2) after the first treatment course.
In the B-mode assessment, we showed that the change in tumor volume following subsequent courses of NAC is not an accurate parameter when analyzed independently; however, in combination with the change in echogenicity, it improved the AUC (0.8 vs. 0.63) and specificity (up to 94%), although with no improvement in sensitivity (64%) after the third NAC course for the non-responding tumors (for RMC ≥ 70%). However, the necrotic areas within the tumor are likely to represent a limitation of this parameter. These areas appear as hypoechoic and they do not affect the change in tumor volume assessed in the B-mode study, which may falsely indicate a lack of response to treatment. However, in studies using MRI in monitoring the response to NAC in BC, Minarikowa et al. showed that, as in the ultrasound evaluation, the appearance of necrosis led to false results [28]. A study showed that higher ADC values before treatment are associated with the presence of necrosis and limited perfusion, and a limited response to treatment should be expected in these tumors as a result. Chu et al. showed that the increase in the ADC parameter value predicts the pCR with a sensitivity of 88% and a specificity of 79% [29]. ADC is a quantitative parameter that measures diffusivity derived from diffusion-weighted imaging (DWI), which is a non-contrast method that evaluates water mobility and tissue cellularity. Increasing ADC values during NAC therapy reflect increased cell lysis and necrosis [30].
Marinovich et al. showed in a meta-analysis that the assessment of the changes in tumors size after NAC assessed by ultrasound and MRI indicates an underestimation of the size of tumors in the ultrasound examination; however, the MRI examination has a tendency to overestimate the size of the residual tumor in relation to histopathological verification [30]. The RECIST criteria used in MRI are not recommended for the assessment of changes in tumor size in ultrasound examination [31,32]. In this study, authors showed that US examination was operator-dependent and characterized by low repeatability. On the other hand, for the next parameter in the B-mode assessment, echogenicity, we showed statistically significant correlation with the RMC value after the third and fourth NAC courses. For the classification of non-responders (RMC ≥ 70%), after the third NAC course, AUC, specificity, and sensitivity were 0.81, 86%, and 64% respectively.
For comparison, in a study on 42 focal lesions, Dobruch-Sobczak et al. analyzed changes in echogenicity after the third NAC course, assuming that tumors did not respond if RMC values were ≥70% and showed sensitivity of 73%; specificity of 87%; PPV of 67%; NPV of 90%; accuracy of 83%; and AUC = 0.69. In the same group, in two parameters analysis using echogenicity and stiffness, an improvement in the statistical parameters was obtained: sensitivity of 82%, specificity of 90%, PPV of 75%, NPV of 93%, accuracy of 88%, and AUC of 0.88. It was noticed that both the stiffness and the hypoechoic nature of the lesions remained unchanged [15].
Therefore, based on the results of statistical analysis, we can assume that changes only in tumor echogenicity, for the group of tumors with RMC ≤ 30%, do not allow the correct assessment of the effects of NAC. In this group, the distinction between residual tumor and fibrosis in the tumor stroma is not possible by analyzing the variability of echogenicity assessed in relation to adipose tissue in the classic B-mode examination. Changes in the residual tumor cells, if present, are variable in pathology examination. More commonly, in microscopic analysis, carcinomas become less cellular and often present as small scattered nests across the tumor bed. After a complete response, only the oedematous vascularized fibroelastotic area of connective tissue with chronic inflammatory cells and macrophages mark the tumor bed [33]. Before treatment, the tumors are rich in malignant cellularity and transmit RF sound better than surrounding breast tissue, and are predominantly hypoechogenic. In contrast, in tumors resistant to NAC, no morphological alteration could be observed, and therefore, the echogenicity was unchanged.
In ultrasound diagnostics, the echogenicity of the B-mode image is a qualitative parameter, even if it is assessed in comparison to the adipose tissue. Although this approach significantly reduces the influence of the ultrasonic scanner settings on the assessment of echogenicity, it remains a subjective assessment by the physician. Image echogenicity is strictly related to the amplitude of the ultrasonic signal scattered in the tissue and the distribution of this amplitude. In the selected area, the assessment of echogenicity is largely dependent on the mean amplitude of the signal. Therefore, the same echogenicity values can be assigned to different amplitude distributions. In our research, we used the KDE (Kernel Density Estimation) method to describe the amplitude distribution of images determined from RF signals, and the KLD (Kullback-Leibler Divergence) method to detect differences between the amplitude distributions of tumor images in subsequent cycles of NAC. The method proposed in this study makes it possible to differentiate the amplitude distributions even when they have the same mean value, which, in relation to the average echogenicity, also enables better tracking of changes taking place in the tissue. Additionally, KLD assessment overcomes the limitations of the subjective assessment of echogenicity based on B-mode images.
Changes in the amplitude distributions in the tumor after particular courses of treatment and the distribution before treatment (KLD 0 ) did not significantly correlate with the RMC for the first and second cycles of NAC, and after the subsequent cycles, the correlation was low (0.3-0.4), but statistically significant. This was reflected during the classification of tumors, which translated into low values of the classification parameters (see Table A3 in Appendix A). Compared to KLD 0 , the KLD 1 parameter correlated better with the RMC.
Applying KLD 1 to predict a poor tumor response (RMC ≥ 70%), the AUC values after the second, third, and fourth NAC cycles were 0.72, 0.82, and 0.81, respectively. The best results were achieved using ∆KLD. After the second course, NAC was able to indicate the non-responding tumors with sensitivity of 84%, specificity of 83%, accuracy of 84%, NPV of 94%, and AUC = 0.85. These values were even higher during the next stages of treatment. The best values after the third treatment course were 82%, 94%, 91%, 94%, and 0.90, respectively. For comparison, the assessment of the effectiveness of treatment based on changes in echogenicity in the B-mode examination allowed us to obtain accuracy, specificity, and NPV exceeding 80% after the third NAC course; however, these values were lower than those obtained for the ∆KLD parameter. This trend was unchanged after the fourth NAC treatment.
In the statistical analysis combining both parameters (∆KLD and ∆Echo) after the second NAC course, all parameters describing the ability to classify remained at the same level as those for classification, which was based only on the ∆KLD parameter. After the third course of NAC, an increase in sensitivity was obtained from 82% to 91%, with a slight decrease in specificity from 94% to 92%, and the values of accuracy and AUC remained at the same high level of 0.91. It is worth noting that this parameter was also a predictor with high potential for the RMC cut-off value of ≤30%. After the second and third courses of NAC, the AUC values obtained were 0.83 and 0.84, respectively. These results suggest that the use of ∆KLD has great potential as a tool for predicting treatment response to NAC.
Note that ∆KLD actually uses two parameters, KLD 0 (n) and KLD 0 (1). Such effective operation of the classification with the use of ∆KLD results from the value of the KLD 0 (1) parameter. The values of this parameter, as shown in the Figure 5, after the first NAC cycle, are higher for tumors with RMC ≥ 70% than for tumors with RMC ≤ 30%. The KLD 0 (1) parameter is in the denominator of the formula describing the ∆KLD parameter.
As a result, dividing by KLD 0 (1) results in an additional elevation of the ∆KLD value for non-responders compared to responders, which improves the classification.
We acknowledge that our study may have some weakness. In this research, postbeamformed RF data were used to avoid the unknown impact of image post-processing algorithms implemented in commercial scanners. RF data are not directly available in common scanners. However, they are always present in the data processing pipeline, and thus QUS methods can be potentially easily implemented on new US platforms.
We only used RMC as a reference standard, which is one of the six residual cancer burdens (RCB) assessed by the pathologist, but in the publication the most common trait referenced by many authors is tumor cellularity.
Another limitation of the study is that only one medical doctor assessed all ultrasound images, however, she has worked for 21 years at the Institute of Oncology and has extensive experience in breast imaging. In B-mode examination, we followed ACR BIRADS-lexicon.

Conclusions
In this study, evaluating the response of BC to NAC after the first four treatment courses, the ∆KLD parameter, which provides quantitative information on changes in image echogenicity, is an accurate predictor of poor response to treatment (RMC ≥ 70%) after the second course of therapy. In the statistical analysis combining ∆KLD and ∆Echo, an increase in sensitivity was obtained without significant deterioration of other statistical parameters.
Our research demonstrates that alterations in tumor echogenicity during NAC treatment are important features in assessing the therapy outcome. The evaluation of changes in echogenicity on the basis of statistical measures applied to the amplitude distributions of the ultrasound image is a particularly promising method. This approach makes it possible to predict the pathological response after NAC. The number of cases used in our research is a limitation in making more general conclusions. Nevertheless, we believe that our results show an important role of quantitative ultrasound in predicting the effects of chemotherapy in breast cancer patients.

Institutional Review Board Statement:
This prospective, single-center study was conducted in accordance with the Declaration of Helsinki, and the protocol was approved by the Ethics Committee of Maria Skłodowska-Curie National Institute of Oncology, Scientific Centre, Warsaw, Poland (Project identification code 49/2018). All participants provided informed consent for inclusion before participation in the study.
Informed Consent Statement: Informed consent was obtained from all subjects involved in the study.
Data Availability Statement: Data supporting reported results can be found in Ultrasound Department, Polish Academy of Science.

Conflicts of Interest:
The authors declare no conflict of interest.      Figure A1. ROC curves for the examined classifiers and four RMC thresholds.