Morphologic complexity of epithelial architecture for predicting invasive breast cancer survival

Background Precise criteria for optimal patient selection for adjuvant chemotherapy remain controversial and include subjective components such as tumour morphometry (pathological grade). There is a need to replace subjective criteria with objective measurements to improve risk assessment and therapeutic decisions. We assessed the prognostic value of fractal dimension (an objective measure of morphologic complexity) for invasive ductal carcinoma of the breast. Methods We applied fractal analysis to pan-cytokeratin stained tissue microarray (TMA) cores derived from 379 patients. Patients were categorized according to low (<1.56, N = 141), intermediate (1.56-1.75, N = 148), and high (>1.75, N = 90) fractal dimension. Cox proportional-hazards regression was used to assess the relationship between disease-specific and overall survival and fractal dimension, tumour size, grade, nodal status, estrogen receptor status, and HER-2/neu status. Results Patients with higher fractal score had significantly lower disease-specific 10-year survival (25.0%, 56.4%, and 69.4% for high, intermediate, and low fractal dimension, respectively, p < 0.001). Overall 10-year survival showed a similar association. Fractal dimension, nodal status, and grade were the only significant (P < 0.05) independent predictors for both disease-specific and overall survival. Among all of the prognosticators, the fractal dimension hazard ratio for disease-specific survival, 2.6 (95% confidence interval (CI) = 1.4,4.8; P = 0.002), was second only to the slightly higher hazard ratio of 3.1 (95% CI = 1.9,5.1; P < 0.001) for nodal status. As for overall survival, fractal dimension had the highest hazard ratio, 2.7 (95% CI = 1.6,4.7); P < 0.001). Split-sample cross-validation analysis suggests these results are generalizable. Conclusion Except for nodal status, morphologic complexity of breast epithelium as measured quantitatively by fractal dimension was more strongly and significantly associated with disease-specific and overall survival than standard prognosticators.


Background
The prognostic assessment of breast cancer is based on factors that determine a patient's relapse risk, and together with predictive factors (e.g., estrogen-receptor status), it is used to make optimal therapeutic decisions regarding adjuvant systemic therapy [1]. Such decisions provide a balance between the potential benefit and associated costs and side effects of treatment [1]. Therefore, it is necessary to have sensitive and specific prognosticators to accurately define risk category for breast cancer.
Currently, the most significant prognosticator for women with breast cancer is axillary lymph node status [1][2][3][4]. For node-positive patients, there is a direct relationship between the number of involved axillary nodes and the risk for distant recurrence [4]. However, despite the usefulness of lymph node status, recommendations for systemic adjuvant chemotherapy are not entirely straightforward. For example, five-year survival rates show that approximately 15% of all node-negative patients with larger tumor sizes (>1 cm) may benefit from systemic adjuvant therapy, but about 85% would survive without it [5]. Furthermore, approximately onethird of node-positive patients are free of recurrence after local-regional therapy [6][7][8].
Other major prognostic risk factors, especially for node-negative patients, are tumor size and histological tumor grade [1][2][3][4]9,10]. For node-negative patients, tumor size is a powerful prognostic factor that is used routinely to make adjuvant treatment decisions [6,11], and tumor grade is primarily used to make decisions for cases in which the tumor sizes are borderline [1,2,5]. Although tumor grade has prognostic value, significant inter-observer variation in grading still exists [12][13][14]. as pathologists are assessing complex histological characteristics in a semi-quantitative manner.
It is known that invasive breast cancer (a malignant neoplasm) demonstrates partial or complete lack of structural organization and functional coordination with surrounding normal tissue [15]. The idea central to this study is that this loss of structural organization and functional coordination manifests itself in the form of an increase in morphologic complexity of the epithelial components at the sub-cellular, cellular, and multi-cellular levels, and the degree of this complexity can be quantified and related to patient outcome. A method that lends itself particularly useful for quantitatively characterizing complex pathological structures at different scales, is based on fractal analysis [16,17]. In this study, we assess the prognostic value of a recently developed novel technique [18] to measure the fractal dimension of segmented histological structures of breast tissue microarray (TMA) cores stained with pan-cytokeratin to highlight the morphology of epithelial architecture.

Patient Characteristics
A total of 408 patients with primary invasive ductal carcinoma (IDC) of the breast were selected retrospectively from the Calgary Regional Hospitals after appropriate ethics approval from the Institutional Review Board (IRB). It should be noted that the IRB did not require patient consent for this study as it was a retrospective study in which many of the patients were deceased and the risk of exposing patient confidentiality was extremely low. Of these, 379 patients had at least one of three TMA cores that was sufficiently stained for fractal analysis. The age range of these patients at diagnosis was 34 to 95 with a mean and median age of 65 and 66, respectfully. Stage information was available for 375 of 379 patients with the following frequency distribution: 225 (60.0%) patients were Stage I, 99 (26.4%) were Stage II, and 51 (13.6%) were Stage III. All patients selected had received adjuvant tamoxifen treatment between 1988 and 2006. Cases were identified with Alberta Cancer Board records of patients who had received tamoxifen treatment without chemotherapy. In summary, the inclusion criterion was any patient who had adequate tissue for TMA construction, and had received adjuvant tamoxifen treatment but no adjuvant chemotherapy.

Sample Preparation
Whole sections stained with Hemotoxylin and Eosin (H&E) were used to select tumor areas for the TMA cores. Fourteen breast TMA blocks containing an average of 94 tissue cores were constructed from formalin-fixed, paraffin-embedded, previously untreated breast cancer tissue. To ensure there was no selection bias, three 0.6 mm cores were chosen randomly from cancerous areas of each donor block to construct the recipient TMA core block, and the Leica RM2235 microtome (Leica Microsystems Inc.) was used to cut 4 μm thick sections from each TMA donor block. In a previous study with prostate cancer specimens, we showed that fractal analyses of specimens stained with pan-cytokeratin provide greater classification performance (benign versus high grade) than serial sections of the same specimens stained with H&E [18]. The reason for this is that pan-cytokeratin isolates and highlights the morphology of epithelial components and excludes structures that do express pathological relevance in the form of morphologic complexity (i.e., connective tissue components). Hence, we stained all the TMA sections with pancytokeratin. This staining was performed using Ventana Benchmark LT. Protease 1 antigen retrieval was used followed by Ventana pre-diluted pan-cytokeratin (cat. No. 760-2135) antibody with an incubation time of 32 minutes. A Ventana ultraview™ DAB detection system was used for detection.

Image Acquisition of TMA Cores
Microscopic images of the TMA cores were acquired with an AxioCam HR digital camera (Carl Zeiss, Inc.) mounted on an optical microscope (Zeiss Axioscope) at a magnification of 10 × objective. The AxioCam HR has pixels of size 6.7 μm × 6.7 μm, which are 1.06 μm × 1.06 μm in apparent size at the combined magnifications of 10 × objective and 0.63 × C-mount optical coupling (optical interface between the microscope and digital camera). The images were taken at the camera's native resolution of 1300 × 1030 pixels, and saved in tagged image file format (tif).

Fractal Analysis to Assess Morphologic Complexity
Unlike our intuitive notion of dimension (i.e., topological dimension), fractal dimension can be a non-integer value, and the greater the morphologic complexity of an object, the higher its fractal dimension relative to its topological dimension ( Figure 1). Fractal dimension quantifies the level of structural complexity by assessing the variation in the level of detail in a structure as the structure is examined at different scales [19]. Hence, it lends itself naturally to characterizing irregular structures that maintain a constant level of complexity over a range of scales.
In this study, we applied an automated fractal analysis technique we developed in previous work [18] to quantify the morphologic complexity of breast epithelium, a pathologically relevant histological feature. In summary, this technique involves the following steps: 1. Application of a histological stain to tissue specimens in order to highlight and isolate the histological structures of interest. In this case, these structures include the outlines of the epithelial components comprising the multi-cellular structures (gland formations), cellular structures (individual cell shapes), and sub-cellular structures (distribution of keratin within the cells and nuclear shape). 2. Image acquisition and background correction of stained specimens. The background correction was done by acquiring a "blank" image (under the same imaging conditions used to acquire the TMA images), and using this "blank" image to subtract the non-uniform background luminance [18]. The resulting background corrected images are converted to grey-scale ( Figure 2). 3. Application of a series of intensity thresholds to convert the grey-scale version of the image specimen into a series of binary images from which histological morphology outlines are derived ( Figure 2). Figure 3 shows a sample magnified region of Figure 2A to illustrate the segmented morphology outlines in more detail. 4. Application of the box counting method [19] (with appropriate spatial scale range -10 to 50 μm) [20] to compute the fractal dimension of each outline image obtained from step 3. 5. Identification of the global maximum from a plot of fractal dimension versus intensity threshold. This maximum corresponds to the fractal dimension of the pathological morphology.
In previous work, we showed that our method of finding the fractal dimension is independent of changes in microscope illumination setting or stain uniformity and intensity [18]. Also, it should be noted that fractal dimension is not affected by magnification as long as the field of view of the specimen image still contains the scale range of the structures of interest over which the fractal dimension was found to be constant.
Our automated fractal analysis method was applied to a total of 1224 TMA cores (3 cores for each of the 408 patient samples). For each patient, the TMA core with the maximum fractal dimension was used for the statistical analysis in this study. The rationale for choosing the maximum fractal dimension from the sampled tissue cores is to reduce the possibility that the other TMA cores from a given patient contain only benign or more highly differentiated tissue. That is, it is expected that the TMA core with the maximum fractal dimension is Figure 1 Both the circle (left) and the Koch snowflake (right) have a topological dimension of 1; however, the fractal dimension (FD) of the Koch snowflake is greater than 1 because it has a more complex morphology than the circle.  representative of the malignant neoplasm that has deviated most from normal cellular/glandular breast morphology, and therefore it is the most probable indicator of abnormal and/or aggressive tumor growth with metastatic potential.
For 379 of the 408 patients (92.9%), fractal dimension was successfully measured in at least one of the three TMA cores generated per patient, and it could not be determined for the remaining 29 patient specimens due to insufficient staining (i.e., less than half of the specimen being stained) or specimen folding. Eight of the 29 patients could not be assessed because all 3 of their TMA cores resulted in a "blank" slide. The breakdown of the number of patients for which the TMA cores were sufficiently stained for fractal analysis was as follows: 36 patients (9.5%) had one evaluable core, 105 patients (27.7%) had two evaluable cores, and 238 patients (62.8%) had three evaluable cores.

Statistical Analyses
For purposes of analyses, it is often useful to convert a measured variable to a categorical variable so as to place patients into graded risk strata. As the particular fractal analysis technique we developed is novel, there are no established cutpoints available. Although several methods exist to determine cutpoints, namely biological determination, data-oriented, and outcome-oriented, there is no single method or criterion to specify which approach is best. For the present analyses, we used a data-oriented approach to select two cutpoints. The first cutpoint was chosen to correspond to the upper quartile (75 th percentile) of the fractal dimension data, and the second cutpoint was chosen as the median of the remaining lower threequarters of the data. Two cutpoints, rather than one, were chosen to assess whether there was a graded relationship between fractal dimension and patient prognosis.
Associations between categorized fractal dimension scores and clinicopathological variables were assessed for statistical significance using a chi-square test. Kaplan-Meier methods were used to estimate 10-year disease-specific and overall survival rates and the logrank test was used to compare the curves for statistical significance. Disease-specific survival was measured from the date of diagnosis to the date of death from cancer or date of last follow-up. Overall survival was measured from the date of diagnosis to the date of death from any cause or date of last follow-up. The above analyses were repeated using Cox proportional hazards regression modeling to assess whether any of the clinicopathological variables influenced the findings. The proportionality assumption was assessed for all covariates using Log-Minus-Log Survival Plots and none violated the assumption. Statistical analyses were performed using SAS 9.2 software (SAS Institute Inc).
The prognostic accuracy of fractal dimension in predicting death from breast cancer and death from any cause was quantified by the area under the curve (AUC) from a receiver operating characteristic (ROC) analysis. Values of AUC range from 0.5 (chance accuracy) to 1.0 (perfect accuracy), with the following intermediate benchmarks: 0.6 (fair), 0.7 (good), 0.8 (excellent), and 0.9 (almost perfect). For the analysis, the predicted probability of outcome from a Cox regression model was considered as a continuum. The actual occurrence of outcome was used as the comparative standard.
A split-sample cross-validation was performed to assess the generalizability of the results [21]. The process consisted of splitting the original sample of 379 patients into a training set of 190 patients and a validation set of 189 patients using random sampling. A regression equation was derived in the training set and the AUC between the observed and predicted response values was calculated. The regression coefficients from the training set were then used to calculate predicted values in the validation set. The AUC between these predicted values and observed values in the validation set was calculated, and is called the cross-validation coefficient. The shrinkage coefficient was calculated as the difference between the AUCs of the training and validation sets. The smaller the shrinkage coefficient, the more confidence one can have in the generalizability of the results. Although there are no clear guidelines regarding the magnitude of shrinkage, except that smaller is better, values less than 0.10 indicate a generalizable model. Given a satisfactory shrinkage coefficient, the data were combined from both sets and a final regression equation was derived based upon the entire sample.
Out of 379 evaluable patients, several had missing data: 15 (9.0%) tumor grades, 4 (1.1%) lymph node status, 15 (4.0%) estrogen-receptor status, and 12 (3.2%) HER-2/ neu status. Rather than excluding these patients from the analyses and reducing the sample size, missing data were imputed using the predicted mean approach in SOLAS 3.0 software (Statistical Solutions, Ltd.). Imputation bias was assessed by re-running all the analyses and excluding any patient with missing data. As the estimates were similar, the results are reported with the imputed data.

Fractal Analysis of the TMA Cores
Fractal dimension scores ranged from 1.08 to 1.97, with a median of 1.62, lower quartile 1.49, and upper quartile 1.75. There was moderate level of relatedness (intraclass correlation = 0.51) among the cores. Using the data-oriented approach to select two cutpoints, fractal dimension values < 1.56 were considered low (N = 141), 1.56-1.75 as intermediate (N = 148), and > 1.75 as high (N = 90). Figure 2 shows representative TMA cores from these fractal dimension categories. One can see from this figure that the classification of TMA cores into low, intermediate, and high fractal dimension categories (A-C) corresponds to the increasing complexity of outline morphology.

Relationship between Fractal Dimension and Standard Prognosticators
The baseline patient characteristics are shown in Table 1. Higher fractal dimension was significantly associated with traditional indicators of poor prognosis, including older age, larger tumour sizes, higher tumour grade, and positive lymph node status. However, fractal dimension was not associated with either estrogenreceptor status or HER-2/neu status.

Fractal Dimension as a Predictor of Outcome
The median patient follow-up was 5.2 years. The 10-year disease-specific and overall survival rates for the entire group of 379 patients were 52.5% and 42.5%, respectively. Patients with higher fractal scores had significantly worse disease-specific survival than those with lower scores (25.0% versus 56.4% versus 69.4%, p < 0.001; Table 2 and Figure 4A). As well, patients with higher scores had significantly worse overall survival (14.2% versus 39.9% versus 67.4%, p < 0.001; Table 2 and Figure 4B). The AUCs for fractal dimension were 0.66 and 0.67 for univariate disease-specific and overall survival, respectively, indicating good levels of prognostic accuracy. As expected, older age, higher grade, and positive lymph node status were significantly predictive of worse outcome, but not the size of the tumour, estrogen-receptor status, or HER-2/neu status (Table 2).

Tumour Grade as a Predictor of Outcome
Tumour grade was derived from the original pathology reports that included between 10 and 30 board-certified cancer pathologists. In contrast to the distinct separation of the disease-specific survival curves for the different fractal dimension categories ( Figure 4A), the disease-specific survival curves for grade 1 and 2 tumours virtually overlaped each other over the entire 10-year follow-up period ( Figure 4C). Also, there is virtual overlap in the overall survival curves of tumour grades 1 and 2 for the first 4-year period ( Figure 4D). These results suggest that tumour grades 1 and 2 do not discriminate patients with respect to 10-year outcome.

Multivariate Analysis
Results from Cox proportional hazards regression showed that fractal dimension remained statistically significant even after adjusting for all clinicopathological variables (Table 3). This result implies that fractal dimension is a strong prognostic factor, even though the multivariate hazard ratio (Table 3) is smaller than the univariate hazard ratio ( Table 2). The AUCs for the 7factor regression models were 0.73 and 0.75 for diseasespecific and overall survival, respectively. These AUCs increased by only 0.07 and 0.08 when six clinical-pathological factors were added to fractal dimension in the multivariate regression model. The small increase in AUCs incidate that the other clinical-pathological factors contribute little to the prognostic accuracy beyond fractal dimension. It is also worth noting that even with the comparison of grades 1 and 2 as one category versus grade 3 tumours, both disease-specific and overall survival were more strongly and significantly associated with fractal dimension than tumour grade.

Split-sample Cross-validation
The generalizability of the aforementioned results was assessed by split-sample cross-validation as described in the statistical analysis section. The results, shown in Table 4 are congruent, not only with each set but also with the results of the entire sample shown in Tables 2  and 3. Specifically, the frequency distribution of low, moderate, and high fractal dimension is similar, as are the 10-year disease-specific and overall survival rates in these three categories. Even with smaller sample sizes, both the training and validation sets still show a pattern of doubling of hazards with higher levels of fractal dimension. The shrinkage coefficients for disease-specific and overall survival were -0.01 and -0.05, respectively, both indicating that fractal dimension is generalizable and that combining data from both sets in the analyses was justified.

Discussion
We previously developed a fractal analysis method to quantitatively measure the morphologic complexity of epithelial architecture [18], and showed a direct association between fractal dimension and breast tumour grade, suggesting that it may be a good surrogate measure of tumour differentiation [22]. In this study we examined the prognostic value of fractal dimension by analyzing 379 specimens from patients with invasive breast cancer, and found that with the exception of nodal status, fractal dimension showed a stronger association with disease-specific survival than standard clinical prognosticators. The potential clinical implications of these results are substantial because to our knowledge, this is the largest and only study of its kind investigating and demonstrating a positive association between the morphologic complexity of breast epithelial architecture (via the fractal dimension metric) and patient outcome. The potential advantages of fractal dimension over conventional tumour grading is that it is a quantitative and reproducible indicator that would be able to provide pathologists with rapid and cost effective high volume analysis from as few as three tissue microarray (TMA) cores per patient. Ideally, a study investigating the value of a potential prognosticator should only involve patients that have not received any form of adjuvant systemic therapy. However, as noted by Mirza et al. [5], such studies are becoming increasingly difficult to perform because systemic therapy is recommended for an ever-widening range of breast cancer patients. Although none of the patients in this study were treated with adjuvant chemotherapy, they were all treated with adjuvant tamoxifen therapy, including the 24 ER-negative patients (note: cases selected for this study where from as far back as 1988 when tamoxifen was occasionally administered to patients with ER-negative tumours). However, even though the patients received a form of adjuvant systemic therapy, the same form of treatment was received by all of the patients leading to the expectation that fractal dimension will be independent of the predictive factor related to tamoxifen therapy (i. e., ER-positive status). Indeed, this appears to be the case, since approximately the same percentage of ERpositive patients are in the low, intermediate, and high fractal dimension groups (Table 1), which likely indicates that tamoxifen therapy has put all of these ERpositive patients on an equal footing. However, another possibility for this result may be that ER status does not affect the morphologic complexity of epithelial architecture. In either case, it may be argued that the use of tamoxifen treated patients in a study investigating the value of a possible prognosticator, although not ideal, does not detract from the ability to assess the prognostic factor's potential relative to other independent prognosticators.
Previous studies have examined the application of fractal analysis for characterizing cancer [23,24] and have shown that fractal dimension can describe the complex pathological structures seen in some cancers; [18,22] however, to our knowledge, our results represent the largest and sole study relating fractal dimension of epithelial architecture to patient outcome. Although we did not use an external patient validation set in this proof of principle study, we employed a data-oriented approach to minimize bias in the selection of cutpoints, as well as, conducting a split-sample cross-validation analysis. This analysis suggests that the results are generalizable, whereby higher fractal dimensions are associated with poorer outcome. This observation demonstrates the high potential of fractal dimension as an image-based prognostic marker, and it is congruent with the notion that malignant breast neoplasms associated with poorer outcome demonstrate partial or complete lack of structural organization and functional coordination with surrounding normal tissue [15]. Furthermore, it implies that changes in the morphologic complexity of architectural components of the neoplasm (i.e., the epithelium) that arise from changes in the

Conclusions
In summary, the results of this retrospective study show that fractal dimension is a promising image analysis marker for the prognosis of IDC of the breast. However, its' prognostic value needs to be confirmed in external validation studies, and ultimately in the context of controlled prospective clinical trials. As a step in this direction, in future work, we will investigate the prognostic value of fractal dimension for defining risk category for Stage I (i.e., lymph node-negative and tumour size ≤ 2 cm in maximum diameter), IDC, ER-positive breast cancer patients that have not received any form of adjuvant systemic therapy. Such a study would be especially valuable because in current clinical practice it is still difficult to identify this subgroup of patients that would benefit most from adjuvant chemotherapy. Also, in future work we will investigate the prognostic and predictive value of combining fractal dimension, a morphological index, with a quantitative analysis of mitotic count, which is a cellular proliferation index that has been shown to be a significant prognostic indicator for node-negative breast cancer [5]. These investigations would provide validation of the significance of morphologic complexity of epithelial architecture in node-negative breast cancer, and explore the possible synergy between morphologic complexity and cellular proliferation. Also, they will bring us closer to the realization of an objective prognosticator that can assist clinicians in making optimal treatment decisions regarding adjuvant systemic therapy for invasive breast cancer.