Geometric Biomedical Computing
Breast cancers can be histologically categorized (graded) based upon their architectural patterns and cellular types. Inaccurate histologic grading can result in inappropriate treatment for a given patient. Computational analysis of breast cancers offers an operator-independent method for histologic grading that should enhance grading reliability. Our approach for automatically and objectively estimating histologic grade is based on image processing and shape analysis of imaged histologic sections. Our work is based on the hypothesis that cellular structures found in breast cancer tumors can be transformed into distinct high-resolution shape distributions using geometric measures from stochastic geometry. The resulting shape distributions define well-populated regions of the associated high-dimensional space. Mapping an unknown breast cancer sample into this high-D space and determining to which region it belongs will allow for the automatic estimation of its histologic grade. See Zhang et al. 2008, Petushi et al. 2011.
In our second project our analysis techniques were used to predict if a breast tumor has metastasized to nearby lymph nodes.
Axillary lymph node metastasis status is still one of the most critical prognostic variables for the breast cancer management decision-making process and patient survival. Metastasis status is determined via examination of a dissected sentinel axillary lymph node. Questions have been raised about the need for routine lymph node dissection. Recently, it has been reported that complete axillary lymph node dissection did not improve survival for patients with small metastatic foci in sentinel lymph nodes. If sentinel lymph node status could be accurately predicted prior to the surgical procedure (lymph node removal), a considerable number of patients with low probability of sentinel lymph node metastases might avoid the procedure altogether, and its associated critical side effects (e.g. swelling, numbness, pain, infection and compromised immunity) and morbidity.
Our approach to predicting breast tumor metastasis is based on image processing and shape analysis of imaged histologic sections at multiple physical scales. We assume that the structure of the nuclear pleomorphisms found in breast cancer tumors can be transformed into distinct high-dimensional shape distributions using geometric measures from stochastic geometry. This geometric space is augmented with information derived from the color variations of the hyperchromatism found in the cancer cell nuclei. We show that shape/color distributions can uniquely capture and characterize the spatial/spectral distribution of neoplastic cells in breast cancer tumors. Once the structure of the cells has been transformed, the dimensionality of the computed distributions is reduced and the resulting feature vectors map into separable regions of the distribution space. Imaging a primary breast carcinoma of unknown status, processing and mapping the image's shapes and colors into feature vectors allow us, via the application of an ensemble-of-classifiers-based approach, to automatically predict a sample's axillary lymph node status, i.e. determine if the cancer has metastasized to nearby lymph nodes by only examining the cellular structures of the primary breast tumor. See Zarella et al. 2015b.
Digital imaging of hematoxylin and eosin (H & E) stained slides has enabled
the application of image processing to support pathology workflows.
However, the intrinsic variability of biological tissue and the vast
differences in tissue preparation protocols lead to significant image
variability that can hamper the effectiveness of these computational tools.
To address this issue we have explored alternative color metrics for image
representation, and transformations to normalize pathology images for
improved standardization and calibration in the anatomic pathology
laboratory.
The
color normalization software tools
developed during our
work (See Zarella et al. 2015a,
Zarella et al. 2017a)
may be accessed
here.
In our latest work we developed a deep-learning-based technique for classification map generation for prostate cancer specimens. See Tan et al. 2021.
Additionally, we have applied an agent-based approach (See Grimes and Breen 2023) to predict metastasis status in breast cancer tissue. See Grimes et al. 2021.
J. Zhang's Research Day 2008 Poster
M.A. Reza's Research Day 2010 Poster
Pathology Informatics 2011 E-Poster
M. Zarella's Discovery Day 2012 Poster
This research is performed in collaboration with the
Johns Hopkins University's Department of Pathology and
Cancer Treatment Centers of America.
Computational pipeline for automated lymph node metastasis prediction based on image analysis of primary breast tumor histology.
Our initial investigations developed image processing techniques that extend the effective penetration depth of OCT imaging into tissue.
This research was performed in collaboration with the
Johns Hopkins University's Department of Pathology and
Cancer Treatment Centers of America.
(A) Single unprocessed OCT frame. (B) Single frame after applying wavelet denoising. (C) Coregistered average of 100 sequential frames.
The GP process is visualized in the following figure.
MPs have also been extended to contain individual coordinate systems. They now detect the orientations of their neighbors and rotate in the direction of the average orientation of nearby MPs while self-assembling.
L. Bai's Research Day 2008 Poster
L. Bai's Research Day 2009 Poster
MPs forming into a diamond shape
MPs forming into an ellipse shape while self-aligning
MPs forming into a diamond shape while self-aligning
MPs forming into a gear shape while self-aligning
Morphogenic Primitives self-organizing into ellipse, diamond, hourglass and gear-like shapes.
Unexpected and interesting shapes and patterns produced during MP evolution.
MPs self-organizing into an ellipse, diamond, hourglass and cross shape, while aligning their orientations.
In the second study, a more detailed investigation was conducted on early stage Drosophila wing discs in order to quantify the change of the geometric parameters of the disc's individual cells during development. Our results showed that the cells' areas constrict, the heights increase and the volumes stay constant during the beginning of the wing disc's formation.
This research was performed in collaboration with the Dr. Frank Jülicher's Biological Physics Department at the Max Planck Institute for the Physics of Complex Systems and Dr. Christian Dahmann's Group at the Dresden University of Technology
3D reconstruction of a Drosophila wing imaginal disc. Individual cells have been color-coded as a function of their apical cross-sectional area. Key units are in microns2.
3D reconstruction of a Drosophila wing imaginal disc. Individual cells have been color-coded as a function of their height. Key units are in microns.
3D reconstruction of a Drosophila wing imaginal disc. Individual cells have been color-coded as a function of their volume. Key units are in microns3.
3D reconstruction of a Drosophila wing imaginal disc at an earlier developmental stage. Individual cells have been color-coded as a function of their apical cross-sectional area.
3D reconstruction of a Drosophila wing imaginal disc at an earlier developmental stage. Individual cells have been color-coded as a function of their height.
3D reconstruction of a Drosophila wing imaginal disc at an earlier developmental stage. Individual cells have been color-coded as a function of their volume.
Plots showing the relationship between wing imaginal disc cell number (developmental stage) and apical face area (in microns2), cell height (in microns) and cell volume (in microns3) for 20 specimens. A line computed with a least squares fit to the data is drawn in each plot. The results show that the cells' areas constrict, the heights increase and the volumes stay constant during the beginning of the wing disc's formation.
One fly behavior that is well established in the field is courtship conditioning. Courting behavior by male flies in Drosophila follows a linear, stereotyped, and well documented set of behaviors, and this behavior is modified by previous sexual experience. After 1 hour of courting a mated female, males suppress their courtship behavior even towards subsequent receptive virgin females for 1-3 hours. This courtship suppression is measured by the Courtship Index (CI), which is calculated by dividing the total amount of time each male fly spends courting by the total duration of a testing period. CI is the standard metric used to assess learning and memory in courtship suppression analysis.
We have developed a computational approach to fly behavior quantification and characterization based on the analysis of videos of courting fruit flies. The approach includes identifying individual flies in the video, quantifying their size (which is correlated with their gender), and tracking their motion, which also involves computing the flies' head directions. Geometric measures are then computed, for example distance between flies, relative orientation, velocities, contact time between the flies, and the time when one fly looks at another. This data is computed for numerous experimental videos and produces high-dimensional feature vectors that represent the behavior of the flies. We formulated a computational equivalent (Computational Courtship Index) of the existing CI, based on the feature vector values, and compared it with CI. Clustering techniques, e.g., k-means clustering, are then applied to the feature vectors in order to computationally group the specimens based on their courtship behavior. Our results show that we are able to reproduce CI values and automatically differentiate between normal and memory/learning defective flies using only the feature vectors derived from our video analysis.
M.A. Reza's Research Day 2011 Poster
This research was performed in collaboration with the Dr. Dan Marenda's Lab and Dr. Aleister Saunder's Lab in Drexel's Biology Department
Steps in the fly identification and tracking process. Calculating a background image, image subtraction, thresholding and filtering to find the flies. Followed by object tracking, that allows us to calculate the geometric quantities from which we derive feature vectors for each video segment.
To date, the GBC Group has investigated three solutions to the contour-based surface reconstruction problem. In the first project, which was led by Dr. Ken Museth of Linköping (Sweden) University, we utilize velocity-adjusted 2D level set contour morphing. With this approach, morphing one contour into the next sweeps out a 3D surface. This is accomplished by equating time in the 2D contour morphing process with the third spatial dimension. A critical aspect of this approach utilizes distance estimates corresponding to the arc lengths of trajectories that connect the adjacent contours in the image plane. These distances, together with a time-of-arrival, are used to estimate the speeds (in contour normal directions) needed to produce a smooth morph when transitioning between sets of contours. Animations of the method in action are here and here.
The second project investigated the effectiveness of Multi-level Partition of Unity (MPU) implicit models to reconstruct surfaces from noisy input contours. Almost all contour-based surface reconstruction techniques exactly interpolate the input contours. MPU implicit surfaces are a type of point set surface that approximates input data within user-specified error bounds. Thus they offer an approach that inherently copes with noisy data in a controllable fashion. The MPU-based reconstruction technique interprets the contour data (pixels in individual images) as points in 3-space. Since MPU implicits also require normal information, it was necessary to develop an algorithm to estimate surface normals from the stacked contours.
A third technique utilizes spline-based 2D distance field interpolation to produce a volumetric representation of the reconstructed surface. Filtering techniques are employed in order to remove the medial axis discontinuities that are found in distance fields. Additionally, monotonicity-constraining natural cubic splines are used to prevent overshoot during interpolation. This third reconstruction approach has proven to be most effective when processing high-complexity, multi-component contours.
J. Marker's Research Day 2006 Poster
This research was performed in collaboration with the Graphics Group of Linköping (Sweden) University, the DUCoM Laboratory for Bioimaging and Anatomical Informatics and the DUCoM Advanced Pathology Imaging Laboratory.
Contour-based surface reconstruction using 2D level set metamorphosis. Pelvis dataset with 35 input contours.
Contour-based surface reconstruction using MPU implicit surfaces. Mouse embryo skin dataset with 186 input contours (46,204 points). Mouse embryo heart dataset with 34 input contours (4,528 points). Mouse embryo stomach dataset with 34 input contours (4,088 points).
Contour-based surface reconstruction using monotonicity-constrained splines to perform 2D distance field interpolation (48 contours).
Contour-based surface reconstruction using monotonicity-constrained splines to perform 2D distance field interpolation. Segmented and classified breast cancer histology image. Color-coded 3D reconstruction of the breast cancer tumor from 9 slices.
The computational model has been extended to simulate the sorting of heterotypic cell populations. The model for studying cell sorting consists of a subset of features from the chemotaxis-based cell aggregation model. The cells in our sorting experiments do not attach to each other, divide or die, but they do age, and emit and respond to chemoattractant gradients. Additionally new features/parameters have been added to each cell type, namely the excretion of a distinct chemoattractant, a chemotactic response rate and a probability of gradient following for each chemoattractant. Our initial in-silico cell sorting experiments only contained two types of cell populations, T1 and T2. Each cell type emits a unique chemoattractant chemical, C1 and C2 respectively. Both types of cells can sense/respond to both types of chemicals and the strengths of these interactions are defined with parameters L1 and L2. A cell's velocity is proportional to the sum of the sensed, adjusted gradients. Another parameter added in the sorting model is P1 and P2, the probability that a cell will follow the gradient of a specific chemoattractant (C1 and C2) during a simulation time step. If a cell does not respond to a gradient it takes a random step. This computational model produces the sorting results included below.
M. Eyiyurekli's Research Day 2007 Poster
M. Eyiyurekli's Research Day 2008 Poster
This research was performed in collaboration with the Drexel Integrated Laboratory for Cellular Tissue Engineering and Regenerative Medicine,
Cell aggregation simulation results.
Cell sorting simulation results.
Cortical bone modeling and drift plays an important role in determining adult bone shape, quantity and quality. With increasing evidence that the bone we grow as children likely affects the health of our bone as we age, it is important to better quantify this process through ontogeny. Moreover, automated means of quantifying these processes can be useful for evaluating the consequences of various health conditions and nutritional deficits on bone growth, important both in modern, orthopedic contexts, and in anthropological contexts to gain insights into the functional adaptations of past populations. In order to investigate this issue the GBC Group is developing a methodology for user-guided segmentation of microCT images of cortical bone cross-sections to (1) discriminate regions of periosteal primary bone apposition that reflect the history of cortical drift in a long bone shaft section and (2) perform measurements that can provide information on bone shape through ontogeny that can then be computationally compared with other the shapes of other bones.
1) A microCT scan of an adolescent bone. 2) Screenshot of program for specifying segmentation parameters. 3) Results of the segmentation. Pores that fall within a specified size range are highlighted in green. 4) Radial lines are used to calculate the thickness of the periosteal primary bone region (between the circular polyline and the bone boundary) as a function of angle.
