Osteosarcoma is the most common primary malignant bone tumor in children and young adults, and its standard treatment follows a neoadjuvant chemotherapy model: patients receive several cycles of chemotherapy before surgical resection, and the resected specimen is then examined by a pathologist to estimate how much of the tumor was killed by the drugs. This histological necrosis percentage has been the cornerstone prognostic indicator in osteosarcoma for four decades, directly informing decisions about whether to continue with the same post-operative chemotherapy regimen or escalate treatment. International series consistently supported necrosis as a predictor of outcome, and the landmark EURAMOS-1 trial enrolled thousands of patients to randomize therapy based on tumor necrosis response.
The measurement problem: Despite its clinical centrality, the histological estimation of tumor necrosis is a manual, time-intensive process. Pathologists must physically section and prepare resected bone tumor specimens, generating as many as 50 histology glass slides per case. The slides are then examined manually, representing only a single plane through a large three-dimensional tumor. The tissue processing protocol for resected osteosarcoma specimens has not been substantially modified since it was first described, and the process is highly dependent on pathologist expertise. Interobserver variability in necrosis estimation is well-documented, with studies reporting meaningful disagreement rates between centers.
Radiological alternatives fall short: Non-invasive imaging approaches, including FDG-PET/CT, diffusion-weighted MRI, and dynamic contrast-enhanced MRI, have been investigated as surrogate biomarkers of chemotherapy response in osteosarcoma, but none has been reliably validated for routine clinical use. Voxel-level interrogation of large three-dimensional tumors is not performed in standard care, where the primary imaging purpose is delineating tumor margins for surgical planning rather than quantifying treatment response.
The critical timing constraint adds to the challenge: histological necrosis estimation is typically available only after 10 weeks of pre-operative chemotherapy, at which point clinically meaningful adjustments to the treatment plan become difficult to implement. This study from UT Southwestern and UT Dallas represents the first fully automated pipeline for classifying osteosarcoma whole slide images into viable tumor, necrotic tumor, and non-tumor regions, directly addressing these long-standing limitations through advances in digital histopathology and automated machine learning.
Archival tumor samples from 50 osteosarcoma patients treated at Children's Medical Center in Dallas between 1995 and 2015 were retrieved from the pathology archives at the time of surgical resection. From these cases, 942 histology glass slides were digitized into whole slide images (WSI) using commercially available scanning technology, with a mean of 19 slides per patient and a range of 4 to 51 slides depending on tumor size and complexity. These WSIs preserve image resolution up to 40X magnification while converting the physical glass slide into a digital file that can be processed computationally.
WSI selection and tile generation: From the 942 digitized images, two expert pathologists selected 40 WSIs specifically chosen to represent the heterogeneity of osteosarcoma morphology and the range of chemotherapy response patterns observed in the full cohort. From each of these 40 WSIs, 30 image tiles of size 1024 x 1024 pixels were randomly generated at 10X magnification, a level chosen in consultation with the pathologists as the optimal balance between spatial resolution and field of view for distinguishing the three tissue classes. Lower magnification would lose cellular detail; higher magnification would introduce irrelevant subcellular information. After excluding tiles falling in non-tissue regions (background), ink-mark regions, and blurry images, 1,144 tiles were retained from an initial set of 1,200.
Pathologist annotation: A custom annotation tool was built in-house specifically for this task, allowing pathologists to review each tile by its unique Tile Identification Number (TIN) and assign one of three class labels: viable tumor (VT), necrotic tumor (NEC), or non-tumor (NT). The 1,144 tiles were divided between the two pathologists, with each tile annotated by a single pathologist. The final distribution reflected the clinical reality of osteosarcoma resection specimens: 536 tiles (47%) were non-tumor, 345 tiles (30%) were viable tumor, and 263 tiles (23%) were necrotic tumor. Annotations were exported to CSV files for subsequent machine-learning training.
Patch generation for deep learning: For the convolutional neural network model, each 1024 x 1024 tile was further subdivided into overlapping 128 x 128 pixel patches, yielding a total of 56,929 image patches. Each patch inherited the class label assigned by the pathologist to its parent tile. This patch-level dataset provided the volume of training examples required for the deep-learning architecture, though it also introduced the assumption that all patches within a tile share the same class label, which may not perfectly hold at tile boundaries.
Before training machine-learning classifiers, the team needed to convert raw image pixels into structured numerical features that capture the tissue properties distinguishing the three classes. They developed two parallel feature generation pipelines: an expert-guided pathway encoding pathologist knowledge about how osteosarcoma tissues look, and an automated pipeline using CellProfiler software to extract a broad set of quantitative image descriptors. The two feature sets are mutually exclusive by design.
Expert-guided features: The eight expert-guided features were derived by translating pathological intuition about what distinguishes viable tumor, necrotic tissue, and normal bone or stroma into programmable image operations. Tiles were converted from RGB to LAB color space, which separates luminance from color information. A threshold of delta >= 64 applied along the a* dimension removed bone, red blood cells, and other false-positive regions. Otsu thresholding then divided the image into foreground and background pixels. Shape segmentation identified nuclei and computed properties including circularity, area, and perimeter. The eight final features covered total nuclei clusters, average clusters per 32x32 window, counts and percentages of foreground vs. background pixels, and average nuclei area and circularity.
CellProfiler features: The automated CellProfiler pipeline generated 53 features through a multi-step process. Color deconvolution using the Unmix Colors module separated each tile into hematoxylin and eosin grayscale components, isolating the nuclear staining from the cytoplasmic background. The hematoxylin image was then analyzed with Haralick texture features computed from a Gray Level Co-occurrence Matrix (GLCM), capturing 13 statistical properties of pixel co-occurrence including angular second moment, variance, sum entropy, sum average, correlation, and Gabor texture response. Otsu thresholding followed by identification of primary objects (nuclei, 30-120 pixel diameter) and secondary objects (cytoplasm, stroma) allowed density-based features to be computed from the Measure Object Neighbors module. Only textural features with high information gain were retained, and features overlapping with the expert-guided set, such as nuclei count, were dropped.
Feature importance comparison: Information gain analysis revealed striking differences between the two feature sets. The top expert-guided feature (total nuclei clusters) had an information gain of 0.2014, while the top CellProfiler features (Texture Gabor, sum entropy, sum variance, variance, sum average, correlation, and angular second moment) all achieved information gain values of 1.8406. Despite this gap, combining both sets outperformed using either set alone. CellProfiler features alone achieved 82.7% accuracy with SVM versus 70.5% for expert-guided features alone, but the combined 61-feature set (53+8) reached 89.9% accuracy, a gain of 7 percentage points over CellProfiler features in isolation. This suggests the expert-guided features encode subtle morphological signals that the automated pipeline misses, even if their individual information gains appear modest.
The combined 61-feature set was used to train and evaluate 13 distinct machine-learning model configurations spanning three algorithmic families: complex decision trees, support vector machines (SVMs), and ensemble learners. This comprehensive evaluation aimed to identify the best-performing model class and configuration for the osteosarcoma classification task, rather than committing to a single approach without empirical comparison. All experiments were conducted within MATLAB's Classification Learner package, which provides standardized implementations of these methods.
Decision trees: Four decision tree configurations were tested, varying the number of splits between 30 and 100, and alternating the splitting criterion between Gini's diversity index and maximum deviance reduction. Decision trees are non-parametric models that partition the feature space into axis-parallel rectangles and assign class labels based on which training examples fall into each region. For each feature, information gain is computed to determine which feature and threshold best separates the three tissue classes at each node. The optimal overall accuracy achieved by complex trees was 80.9%.
Support vector machines: Six SVM configurations were evaluated, varying the kernel function among quadratic, cubic, and Gaussian radial basis function kernels, and alternating between one-vs-one and one-vs-all multi-class decomposition strategies. SVM models identify a hyperplane in the feature space that maximally separates training examples from different classes, with slack variables allowing some misclassification in the soft-margin formulation. The best-performing SVM used a cubic polynomial kernel (degree 3, designated SVM3). With 5-fold cross-validation on the combined feature set, SVM3 achieved the highest overall accuracy of 89.9%, with class-specific accuracies of 91% for viable tumor, 87% for necrotic tumor, and 91% for non-tumor.
Ensemble learners: Three ensemble configurations were tested using bagging (bootstrap aggregating) and boosting (AdaBoost) on decision tree base classifiers. Bagging reduces variance by training on multiple bootstrap samples and averaging predictions. Boosting reduces bias by iteratively reweighting training examples, placing more emphasis on those misclassified in prior rounds. Ensemble learners achieved overall accuracies ranging from 81.1% to 86.8%, placing them between decision trees and SVMs. The fact that SVMs consistently outperformed both tree-based approaches on this dataset reflects the advantage of margin-based classification in high-dimensional feature spaces with a relatively small number of training examples.
In parallel with the machine-learning pipeline, the team built a custom convolutional neural network (CNN) architecture inspired by AlexNet and LeNet, designed to operate on the 128x128 pixel patches generated from each annotated tile. Unlike the machine-learning models, which required explicit hand-crafted features as input, the CNN learns feature representations directly from raw pixel values through the training process. This end-to-end approach eliminates the need for a separate feature engineering stage, potentially capturing morphological patterns that do not map cleanly onto the feature categories captured by the CellProfiler or expert-guided pipelines.
Network architecture: The CNN consists of three alternating convolution and max-pooling layer pairs, followed by two fully connected multilayer perceptrons (MLPs). The input layer accepts raw RGB pixel values from 128x128 pixel patches. Each convolution layer computes dot products between input values and learned filter weights, extracting local features at increasing levels of abstraction. The first two convolution-pooling pairs use 5x5 filter kernels with a downsampling scale of 2; the third pair uses 3x3 filters. After the three convolution-pooling stages, the spatial dimension is reduced to 62x62x4 before entering the fully connected layers. The final MLP layer has three output neurons producing probability values for each of the three tissue classes (VT, NEC, NT).
Training configuration: The network was trained using stochastic gradient descent with Nesterov momentum (gamma = 0.9). All weights were initialized from a Gaussian distribution N(0, 0.01). Hyperparameters were set by the "babysitting approach," meaning manual tuning across multiple training runs guided by validation loss convergence. The final configuration trained for 20 epochs with a learning rate of 10^-3 and a batch size of 100. Two data augmentation strategies were applied to reduce overfitting and increase effective training set size: geometric augmentation (rotation by 90, 180, and 270 degrees, plus horizontal and vertical flips, yielding 8 transformations per patch) and PCA-based color perturbation that alters the principal color components of each patch while preserving structural features.
Patch preprocessing: Before training, patches were filtered using the Euler value criterion to remove those consisting primarily of background (white) pixels. An Euler threshold of e >= -3, determined by testing values from -15 to 15 on a 1,000-patch validation set containing 300 background-only examples, was applied to select only informationally significant patches. Remaining patches were normalized from RGB to LAB color space with normalization of the L (lightness) channel to reduce variability introduced by differences in slide staining intensity and scanner illumination across the 942-slide dataset spanning 20 years of archival material.
Both the best machine-learning model (SVM3 with cubic kernel) and the CNN were evaluated on the same held-out test set using an 80-20 stratified split, preserving the original class proportions (47% non-tumor, 30% viable tumor, 23% necrotic tumor) in both training and test partitions. The 80-20 split yielded 914 training tiles and 230 test tiles for the tile-level evaluation, with the CNN additionally evaluated at the patch level on its 56,929-patch dataset split using the same 80-20 ratio.
Machine-learning performance: SVM3 achieved an overall accuracy of 89.9% on 5-fold cross-validation and maintained strong class-specific performance: 91% for viable tumor, 87% for necrotic tumor, and 91% for non-tumor. Necrotic tumor was the most challenging class across all model families, consistently yielding lower class-specific accuracy than viable tumor or non-tumor, reflecting the morphological heterogeneity of necrotic tissue where cell death patterns vary widely depending on the tumor type and duration of response.
Deep-learning performance: The CNN achieved an overall patch-level accuracy of 93.3%, with class-specific patch accuracies of 95.3% for viable tumor, 92.7% for necrotic tumor, and 91.9% for non-tumor. At the tile level (aggregating patch predictions back to the 1024x1024 parent tile), tile-level accuracies were 92.6% for viable tumor, 91.5% for necrotic tumor, and 89.5% for non-tumor. The CNN outperformed the SVM on both overall accuracy and on the most clinically critical class distinction, viable versus necrotic tumor.
ROC and VUS analysis: The hierarchical receiver-operator characteristic was computed in two stages: first, discriminating tumor from non-tumor, then conditionally discriminating necrotic from viable tumor. For the binary tumor vs. non-tumor task, both models achieved AUC values of 0.98 or higher, indicating near-perfect separation of tumor from background tissue. For the three-class problem, the Volume Under the ROC Surface (VUS) was computed as the multi-class generalization of AUC. Perfect classification corresponds to VUS = 1.0, while uninformed random classification corresponds to VUS = 0.167. The CNN achieved VUS = 0.959, outperforming the SVM at VUS = 0.922. Both values indicate highly discriminative models, and their agreement across two fundamentally different algorithmic approaches provides strong internal validation of the pipeline's robustness.
Beyond tile-level accuracy metrics, the clinical value of the pipeline depends on how well the trained models can be applied to analyze complete whole slide images, producing outputs that pathologists can interpret and use to estimate the percent tumor necrosis. The authors implemented a visualization step in which trained models are applied to systematically generated tiles spanning an entire WSI, and the classification outputs are mapped back to the spatial location of each tile within the original slide to produce a color-coded tumor-prediction map.
Tile generation from test WSIs: For each test WSI, 1,000 to 1,800 image tiles of size 1024x1024 were generated to provide comprehensive coverage of the tissue sections. Each tile was independently classified by both the SVM3 and CNN models, with the predicted class label (viable tumor, necrotic tumor, or non-tumor) recorded for each tile's spatial position. This systematic tiling approach ensures that no tissue region is missed, unlike manual pathological review which samples representative areas.
Color-coded visualization: The classification output for each tile was mapped back to the original WSI at eye-fit level, creating a color-coded overlay where viable tumor regions are displayed in red, necrotic tumor regions in green, and non-tumor regions in blue. The resulting tumor-prediction map provides a bird's-eye view of the spatial distribution of viable and necrotic tissue across the entire slide, making it straightforward to visually estimate the proportion of necrosis and identify areas where the two tissue types are intermingled, which is common in osteosarcoma with heterogeneous response.
Clinical pathway to necrosis percentage: The tumor-prediction map is designed as the penultimate step in a complete tumor assessment pipeline. The final step, not yet implemented in this study but explicitly planned as future work, would automatically compute the percentage of tumor necrosis by calculating the ratio of necrotic tumor tiles to total tumor tiles (necrotic plus viable) across the map. This automated necrosis percentage, combined with correlation to co-registered MRI sequences, represents the complete vision of the pipeline from digitized histology slide to clinical prognostic biomarker.
Interobserver variability in annotation: The training labels assigned to each tile were determined by a single pathologist per tile, with no second-opinion confirmation or consensus annotation process. The authors acknowledge that pathological evaluation of tissue samples is prone to interobserver variability, and that some of the expert-guided features depended directly on pathologist-identified characteristics. This variability in the ground truth labels propagates into the trained models: if two pathologists would disagree on the class of a tile, the model can only learn one of those interpretations. The parallel CNN approach was partly developed to address this, since it functions independently of explicit pathologist-guided feature definitions and learns directly from image pixels.
Dataset size and single-center origin: While 50 patients and 942 digitized slides represents a substantial archival dataset, the 40 WSIs selected for annotation represent a relatively limited training set for deep-learning by contemporary standards. All cases originated from a single institution (Children's Medical Center, Dallas), meaning the trained models have not been exposed to the staining protocol variability, scanner differences, and patient population diversity encountered across the many centers where osteosarcoma patients receive care. External validation at independent sites with their own pathology laboratories and slide preparation protocols is the essential next step before any clinical deployment.
Tile-level annotation assumption: Each 1024x1024 tile received a single class label, which is then inherited by all 128x128 patches derived from it. In practice, tiles may contain mixed regions spanning two or even all three tissue classes, particularly at boundaries between viable and necrotic tumor or between tumor and non-tumor stroma. Assigning a single label to such mixed tiles introduces label noise that reduces classification performance and may cause the models to learn from ambiguous examples. Pixel-level or sub-tile annotation, while far more labor-intensive, would produce cleaner training labels and likely improve model performance on boundary cases.
Future directions: The authors outline a three-part future agenda. First, the pipeline should be extended to automatically compute the percentage of tumor necrosis from the tumor-prediction map, enabling a direct comparison with manual pathologist estimates. Second, digital histology data should be co-registered with enhanced MRI sequences (diffusion-weighted and dynamic contrast-enhanced) at matching anatomical planes, creating a multimodal dataset that can test whether combining imaging and histological AI predictions outperforms either modality alone. Third, the framework should be adapted for other tumor types, since the modular architecture of the pipeline, from WSI tiling through feature extraction through classification to visualization, is not osteosarcoma-specific and has been demonstrated to generalize in the machine-learning literature for other carcinomas.