Defining a Radiomic Response Phenotype: A Pilot Study using targeted therapy in NSCLC

Medical imaging plays a fundamental role in oncology and drug development, by providing a non-invasive method to visualize tumor phenotype. Radiomics can quantify this phenotype comprehensively by applying image-characterization algorithms, and may provide important information beyond tumor size or burden. In this study, we investigated if radiomics can identify a gefitinib response-phenotype, studying high-resolution computed-tomography (CT) imaging of forty-seven patients with early-stage non-small cell lung cancer before and after three weeks of therapy. On the baseline-scan, radiomic-feature Laws-Energy was significantly predictive for EGFR-mutation status (AUC = 0.67, p = 0.03), while volume (AUC = 0.59, p = 0.27) and diameter (AUC = 0.56, p = 0.46) were not. Although no features were predictive on the post-treatment scan (p > 0.08), the change in features between the two scans was strongly predictive (significant feature AUC-range = 0.74–0.91). A technical validation revealed that the associated features were also highly stable for test-retest (mean ± std: ICC = 0.96 ± 0.06). This pilot study shows that radiomic data before treatment is able to predict mutation status and associated gefitinib response non-invasively, demonstrating the potential of radiomics-based phenotyping to improve the stratification and response assessment between tyrosine kinase inhibitors (TKIs) sensitive and resistant patient populations.


1
Uni / 1D Maximal diameter of tumor. Calculated by multiplying the longest line inside the tumor on an axial plane by the image resolution in x or y direction. 2 Volume / 3D Volume of tumor. Calculated by multiplying the number of tumor voxels by the image resolutions in x-, y-and zdirections.

3-12
First Order Statistics / 2D and 3D Features in this group are derived from the intensity histogram of tumor. The following 5 well-known first order statistics are calculated in both 3D and 2D, totaling 10 features in this group.
• Mean: average intensity of tumor • SD: standard deviation of tumor intensities • Skewness: a measure of intensity symmetry • Kurtosis: a measure of flatness of tumor intensity relative to a normal distribution • Peak Position: the peak position of the histogram [13][14][15][16][17][18][19][20][21] Shape / 2D and 3D This group of features describes the shape properties of the tumor. In 2D, the following three shape features are extracted: • Roundness Factor (RF): a measure of circularity of a tumor's profile on the 2D image. It is defined as, where Area denote the area of tumor and Perimeter the length of tumor contour on the 2D image. • Eccentricity: a measure specifying how close an ellipse (i.e., a tumor's profile on the 2D image) is to a circle. It is defined as: Eccentricity = c / a where, c is the distance from the center to a focus and a is the distance from that focus to a vertex. Shape index is used to intuitively characterize the local shapes of tumor surface [1]. The shape indexes, SI(1) -SI (9), are numbers between [-1, 1] and respectively imply the amount of the following 9 shapes on tumor surface: 1]: spherical cap Since the sum of the 9 scaled shape indexes equals to 1, SI(1) is excluded from the analysis to reduce the redundancy. This group has 8 features.

30-35
Sigmoid Function / 3D To quantify tumor margins, Sigmoid Function is used to fit density change along a sampling line drawn orthogonal to the tumor surface [2]. Each sampling line, going through one voxel on the tumor surface, has a certain length (5mm in this work) inside and outside the tumor. The Sigmoid Function is defined as, where the fitting parameter A, B and C respectively specify the amplitude, slope and offset of the curve. Wavelets / 2.5D In this study, two types of Wavelets transform are implemented. One is discrete wavelet transform (DWT) [3] and the other is discrete stationary wavelet transform (SWT) [4]. For each transform, Wavelets images/features are calculated at three decomposition levels, resulting in 18 features in this group. DWT: Taking a M x N image I(m, n) as an example, the first level DWT decomposition can be briefly described as the following. First, a low-pass and a high-pass filter ('Coiflets1' wavelet filter) are applied to the original image vertically followed with a vertical down-sampling by a factor of 2. Then the two filters are applied to the processed image horizontally followed by a horizontal down-sampling by a factor of 2. This results in 4 subimages that are known as the low-pass approximation horizontal detail H(m, n) and diagonal detail D(m, n). The second (and (third) level DWT decomposition repeats the above procedure but with the average image generated at the first (and second) level decomposition.
In this study, 9 wavelet features are defined as the Energy of each detailed sub-images. Let N i be the number of pixels of a sub-image at level i (i=1, 2).

At the first DWT decomposition level,
• DWT-H:

At the second DWT decomposition level,
• DWT-LH: • DWT-LH: • DWT-LV: To overcome the translation-variance in DWT, SWT performs up-sampling on wavelet filters instead of downsampling on the image. Thus, at each scale of SWT decomposition, the image is convolved by up-sampled filter coefficients and remains the original size. There are seven corresponding SWT features (detail not included).

54-62
Edge Frequency / 2.5D Edge Frequency features, obtained from images processed by an edge operator (in this work, it is a 2D Robert's edge operator), characterize variation of the intensity gradient inside a tumor [5]. The 2D Robert's edge operator is defined as follows: denotes a tumor pixel intensity at the location (i, j), and is the distance between the pixel ( , ) and its neighboring pixel. The three features, Mean, Coarseness and Contrast, are computed from the processed gradient images at the distances of one, four and nine pixels, which totals 9 features. Equations to define the Coarseness and Contrast are the same as the ones defined in the GTDM feature group, with S(i) replaced by gradient(d).

63-64
Fractal Dimension / 2.5D Fractal Dimension provides a statistical index to characterize the complexity of an image [6]. Basically, the Fractal Dimension describes the relationship between the changes in a measuring scale and the measurement results at the scale. In this work, a 3D box-counting algorithm [7] is adopted to calculate the Fractal Dimension to quantify tumor intensity homogeneity. Fractal Dimension Mean and Standard Deviation are computed.

65-79
Gray-Tone Difference Matrix (GTDM) / 2.5D Neighborhood GDTM features are defined based on gray-tone (i.e., image intensity) difference between a pixel and its neighborhood [8]. Let ( , ) be an image pixel that has the gray-tone of and is located at ( , ).
The average gray-tone over a neighborhood centered at, but excluding ( , ), is otherwise where ! is the set of all pixels having the gray-tone of . For an × image, let ! be the probability of occurrence of gray-tone value i, ! be the highest graytone value in the image and ! be the total number of different gray-tone values in the image. The GDTM features are defined as, • Coarseness: In this study, ! = 256. The distance of the neighboring pixels is 1.

80-99
Gabor / 2.5D Gabor filters are linear filters designed for detecting edges at different directions and width [9]. It is an oriented Gaussian function modulated by a sinusoidal wave. The Gabor Energy feature is defined as the sum of the square of intensity over all tumor pixels on the images processed by Gabor filter. In this work, the Energy is calculated at 4 wavelengths (w=3, 5, 7 and 9 pixels) and Laplacian of Gaussian (LoG) is a combined filter, i.e., a Gaussian smoothing filter followed by Laplacian, a differential operator [11]. The definition of a 2D LoG is: In this study, the following three LoG features of Mean, Uniformity and Entropy are calculated from the LoG filtered (processed) image, ( ), at four σ levels σ ϵ 0, 0.5, 1.5, 2.5 . σ = 0 (s1; no smoothing) and σ = 2.5 (s4). This totals 12 LoG features in this group.
• LoG Mean Gray Intensity (MGI): where N is the number of object pixels, P(i) is the probability of pixels with a gray-level of in the LoG pre-processed image, and is the maximal value of the pre-processed image.

126-130 Run-length / 3D
The Run-length features are developed to characterize tumor heterogeneity by counting the number of maximum contiguous pixels / voxels having an identical intensity along a line [12]. Let , be the number of primitives of all directions having length of and intensity of , the volume of tumor, ! the maximum run-length, and = 256 the number of image intensity bins. The total number of run-lengths is therefore, In this work, the following five Run-length features are used.
• Run_SPE: Short primitives emphasis = Spatial Correlation features assess linear spatial relationships between texture primitives [13]. Let , be an image pixel's intensity at the location (x, y) in a tumor, the distance between two pixels, the area of the tumor, ! the area of the tumor after shrinking with a distance of pixels. Then, In this study, the spatial correlations are calculated at the distances of =1 and d=4 pixels. 133-183 Gray-Level Cooccurrence Matrix (GLCM) /3D This feature class quantifies textures by creating a new matrix, called GLCM, which is based on the frequency of image pixel pairs possessing particular intensity values at a certain direction and distance [14]. In this work, GLCMs are generated at 13 directions and three distances (d=1, 4, 9 voxels). At each distance, the final GLCM is the average of the 13 GLCMs. For each of the three GLCMs, the following 17 standard statistical features are derived to characterize tumor's homogeneity, contrast, entropy, etc.