# Comparative Multifractal Analysis of Dynamic Infrared Thermograms and X-Ray Mammograms Enlightens Changes in the Environment of Malignant Tumors

^{1}Laboratory of Physical Foundation of Strength, Institute of Continuous Media Mechanics UB RAS, Perm, Russia^{2}CompuMAINE Laboratory, Department of Mathematics and Statistics, University of Maine, Orono, ME, USA^{3}Université Lyon, Ecole Normale Supérieure de Lyon, Université Claude Bernard Lyon 1, Centre National de la Recherche Scientifique, Laboratoire de Physique, Lyon, France^{4}Laboratoire Ondes et Matière d'Aquitaine, Centre National de la Recherche Scientifique, Université de Bordeaux, UMR 5798, Talence, France^{5}Department of Therapeutic and Propedeutic Dentistry, Perm State Medical University, Perm, Russia

There is growing evidence that the microenvironment surrounding a tumor plays a special role in cancer development and cancer therapeutic resistance. Tumors arise from the dysregulation and alteration of both the malignant cells and their environment. By providing tumor-repressing signals, the microenvironment can impose and sustain normal tissue architecture. Once tissue homeostasis is lost, the altered microenvironment can create a niche favoring the tumorigenic transformation process. A major challenge in early breast cancer diagnosis is thus to show that these physiological and architectural alterations can be detected with currently used screening techniques. In a recent study, we used a 1D wavelet-based multi-scale method to analyze breast skin temperature temporal fluctuations collected with an IR thermography camera in patients with breast cancer. This study reveals that the multifractal complexity of temperature fluctuations superimposed on cardiogenic and vasomotor perfusion oscillations observed in healthy breasts is lost in malignant tumor foci in cancerous breasts. Here we use a 2D wavelet-based multifractal method to analyze the spatial fluctuations of breast density in the X-ray mammograms of the same panel of patients. As compared to the long-range correlations and anti-correlations in roughness fluctuations, respectively observed in dense and fatty breast areas, some significant change in the nature of breast density fluctuations with some clear loss of correlations is detected in the neighborhood of malignant tumors. This attests to some architectural disorganization that may deeply affect heat transfer and related thermomechanics in breast tissues, corroborating the change to homogeneous monofractal temperature fluctuations recorded in cancerous breasts with the IR camera. These results open new perspectives in computer-aided methods to assist in early breast cancer diagnosis.

## 1. Introduction

The past 30 years have seen the emergence in cancer biology of the concepts of local microenvironments and stem cell niches not only as key regulators of tissue specificity, homeostasis maintenance and tumor control, but also as active players in cancer initiation and progression to metastasis (Bissell et al., 2005; Bissell and Hines, 2011; Maguer-Satta, 2011; Lu et al., 2012). Important progress was made in the understanding of the permanent dialogue that exists between stem cells and their microenvironment (Fuchs et al., 2004; Li and Li, 2006; Moore and Lemischka, 2006; Morrison and Spradling, 2008). The role of the niche in normal tissues is critical in regulating asymmetric cell division and stem cell fate (Fuchs et al., 2004; Morrison and Kimble, 2006; Ho and Wagner, 2007; Marthiens et al., 2010). On the one hand, it participates to maintain stem cell quiescence until self-renewal divisions required to sustain the stem cell pool from exhaustion and alteration, and, when needed, to facilitate the proliferation necessary to respond to some physiological demand or injury (Morrison and Kimble, 2006; Morrison and Spradling, 2008; Trumpp et al., 2010). On the other hand, it contributes to drive cell differentiation processes during lineage specification and organ development (Ho and Wagner, 2007; Marthiens et al., 2010; Trumpp et al., 2010). Upon asymmetric division, cell commitment is distributed upon two daughter cells; one of them keeps stem cell properties while the other one is driven toward a more differentiated stage to respond, match, and adapt to the surrounding tissue (Marthiens et al., 2010). But an imbalance between stem cell activation and differentiation can lead to the generation of damaged stem cells by either overexpanding the stem cell pool or a failure in stem cell differentiation (Demicheli, 2001; Bissell et al., 2005; Bissell and Labarge, 2005; Flynn and Kaufman, 2007). In homeostatic conditions, the stromal environment of the niche, including fibroblast, vasculature and immune cells as well as interstitial extracellular matrix (ECM), can biochemically and biomechanically control the so-called cancer stem cells to maintain tissue architecture and integrity (Bissell and Hines, 2011; Maguer-Satta, 2011; Lu et al., 2012). This explains why many occult tumors can lie dormant or evolve very slowly (Demicheli, 2001; Bissell et al., 2005; Bissell and Labarge, 2005; Faraldo et al., 2005; Li and Neaves, 2006; Moore and Lemischka, 2006; Flynn and Kaufman, 2007; Tysnes and Bjerkvig, 2007). There is increasing evidence that the destabilization of tissue homeostasis originates from the alteration of the tumor microenvironment that not only affects the behavior of cancer stem cells, stroma and various epithelial cells, but also contributes to transform the niche into a tumorigenic microenvironment that further facilitates the process of oncogenic transformation, tissue invasion and metastasis evasion during cancer progression (Bissell et al., 2005; Lee and Herlyn, 2007; Rønnov-Jessen and Bissell, 2009; Bissell and Hines, 2011; Maguer-Satta, 2011; Lu et al., 2012). By favoring the survival and proliferation of cancer stem cells, this altered cancer promoting environment can maintain the quiescence of these cells and make them resistant to treatment, creating the possibility of regenerating a tumor (Demicheli, 2001; Hall et al., 2007; Besançon et al., 2009; Bissell and Hines, 2011; Maguer-Satta, 2011). The ongoing progress in understanding the fundamental role of cancer stem cells within their niche raises very challenging issues in cancer assessment, diagnosis and therapy. There are emerging strategies and recent clinical trials relevant to microenvironmental therapies (Hall et al., 2007; Besançon et al., 2009; Bissell and Hines, 2011; Maguer-Satta, 2011). In this work, we combine the analysis of dynamic IR thermograms and X-ray mammograms to show that changes in the environment of a malignant breast tumor can be detected with commonly used non invasive screening techniques.

Breast cancer is one of the major causes of death among women worldwide (Siegel et al., 2015). Clinical studies have demonstrated that survival is significantly improved if the breast anomalies are detected as early as possible (Lee, 2002). Despite some criticism of the use of screening techniques due to overdiagnosis (Jørgensen and Gøtzsche, 2009; Fenton et al., 2011), early detection remains the best strategy for improving prognosis and providing less invasive options for both specific diagnosis and treatment (Ganesan et al., 2013; Jalalian et al., 2013). In a preliminary study (Gerasimova et al., 2013, 2014), we showed that skin temperature dynamics recorded with an infrared (IR) camera displayed qualitative changes around malignant breast tumors. Using a wavelet-based time-frequency method, we demonstrated that temperature temporal fluctuations superimposed on the cardiogenic and vasomotor perfusion oscillations are not instrumental noise, but contain physiological information that can be exploited to anticipate the transition to malignancy. The observed drastic simplification from multifractal (continuous change of statistics across time-scales) to homogeneous monofractal (invariant statistics across time-scales) skin temperature fluctuations in malignant tumor foci, was hypothesized to be the signature of blood vessels and other tissues showing signs of aberration and architectural disintegration, likely affecting heat transfer and related thermomechanics inside the breast. The aim of the present study is to strengthen this interpretation by showing that this disorganization of the breast tissue architecture and vascular network in the environment of a malignant tumor can also efficiently be detected with X-ray mammography, the golden standard for breast cancer screening detection (Nass et al., 2001; Bronzino, 2006). Using a wavelet-based space-scale method (Arneodo et al., 2003), we show that the long-range correlations and anti-correlations respectively observed in the roughness fluctuations of the mammograms of dense and fatty normal breasts, strikingly vanish in the malignant tumor region. Note that most existing computer-aided diagnostic (CAD) methods (Karahaliou et al., 2008; Ayer et al., 2010; Tsai et al., 2011; Häberle et al., 2012; Meselhy Eltoukhy et al., 2012) are designed for texture analysis or feature extraction with the prerequisite that the background roughness fluctuations of normal breast mammograms are homogeneous and uncorrelated. In contrast, we propose characterizing correlations in mammogram roughness fluctuations via the estimate of a density fluctuation index *H* as an innovative and effective discrimination method that will assist in early breast cancer detection.

## 2. Methods of Analysis

The wavelet transform (WT) is a mathematical microscope (Arneodo et al., 1988, 1995, 2008; Muzy et al., 1991, 1994) suitable to the analysis of complex non-stationary time-series, such as those found in genomics (Nicolay et al., 2007; Arneodo et al., 2011; Audit et al., 2013) and physiological systems (Ivanov et al., 1999, 2001; Goldberger et al., 2002; Richard et al., 2015), thanks to its ability to be blind to low-frequency trends in the analyzed signal Σ(*t*). Its generalization in two dimensions (2D) is used in image processing (Mallat, 1998; Arneodo et al., 2003, 2008; Antoine et al., 2008), including applications in cellular biology (Khalil et al., 2007; Snow et al., 2008; Goody et al., 2010; Grant et al., 2010) and medicine (Kestener et al., 2001; Khalil et al., 2009; Batchelder et al., 2014).

### 2.1. The 1D Wavelet Transform (Time-Frequency Analysis)

The WT is a time-scale decomposition method which consists in expanding signals in terms of wavelets constructed from a single function, the “analyzing wavelet” ψ, by means of translations and dilations. The WT of a real-valued function Σ(*t*) is defined as (Mallat, 1998)

where *t*_{0} is a time parameter and *a* (> 0) a scale parameter (inverse of frequency). By choosing a wavelet ψ whose *n* + 1 first moments are zero [∫*t*^{m}ψ(*t*)*dt* = 0, 0 ≤ *m* ≤ *n*], one makes the WT microscope blind to order-*n* polynomial behavior, a prerequisite for multifractal fluctuations analysis (Muzy et al., 1991, 1994; Arneodo et al., 1995, 2008).

### 2.2. The 2D Wavelet Transform (Space-Scale Analysis)

With an adapted analyzing wavelet, one can reformulate Canny's multiscale edge detection in terms of a 2D wavelet transform (Mallat, 1998). The underlying strategy consists in smoothing the discrete image data by convolving it with a filter prior to computing the gradient of the smoothed image. Let us define two wavelets as the partial derivatives with respect to *x* and *y* of a 2D-smoothing function ϕ(*x*, *y*) (Mallat, 1998; Arneodo et al., 2003):

The WT of *I*(*x*, *y*) with respect to ψ_{1} and ψ_{2} has two components and therefore can be expressed in a vectorial form:

where **x _{0}** is a space parameter and

*a*(> 0) a scale parameter. From Equation (3), we can compute the modulus and argument of the WT:

and

### 2.3. The Wavelet Transform Modulus Maxima Method

The wavelet transform modulus maxima (WTMM) method was originally developed to generalize box-counting techniques (Arneodo et al., 1987) and to circumvent the limitations of the structure function method to perform multifractal analysis of 1D velocity signal in fully-developed turbulence (Muzy et al., 1991, 1993, 1994; Bacry et al., 1993; Arneodo et al., 1995). It proved to be a reliable method to estimate scaling exponents and multifractal spectra (Muzy et al., 1994; Delour et al., 2001; Audit et al., 2002). The WTMM method was generalized in 2D for the multifractal analysis of rough surfaces (Arneodo et al., 2000, 2003; Decoster et al., 2000; Roux et al., 2000), and for the analysis of 3D scalar and vector fields (Kestener and Arneodo, 2003, 2004; Arneodo et al., 2008). Successful applications have already been reported in various area of fundamental research (Muzy et al., 1994; Arneodo et al., 1995, 2002, 2008, 2011). In the context of present study, the 1D WTMM method demonstrated the multifractality in physiologic dynamics and its breakdown with disease (Ivanov et al., 1999, 2001; Goldberger et al., 2002), whereas the 2D WTMM method was used to detect microcalcifications and shown to have great potential to assist in cancer diagnosis from digitized mammograms (Kestener et al., 2001; Arneodo et al., 2003; Batchelder et al., 2014).

In 1D, the WTMM method consists in computing the WT skeleton defined, at each fixed scale *a*, by the local maxima ${L}$(*a*) of the WT modulus ${M}$(*t*, *a*) = |*W*_{ψ}(*t*, *a*)|. These WTMM are positioned across scales on curves called maxima lines (Figure S1). In 2D, these WTMM lie, for a given scale, on connected chains called maxima chains (Figures S2A–C). Considering the points along these maxima chains where ${M}$(**x**, *a*) is locally maximum, we define the so-called WTMMM. The WTMMM are then linked through scales to form the WT skeleton (Figure S2D). Along these maxima lines the WTMM (resp. WTMMM) behave as *a*^{h(t)} (resp. *a*^{h(x)}), where *h*(*t*) (resp. *h*(**x**)) is the Hölder exponent characterizing the singularity of Σ (resp. *I*) at time *t* (resp. position **x**) (Muzy et al., 1994; Arneodo et al., 2003, 2008). The multifractal formalism amounts to quantify statistically the contributions of each Hölder exponent value via the computation of the singularity spectrum defined as the fractal dimension *D*(*h*) of the set of points *t* (resp. **x**) where *h*(*t*)(resp. *h*(**x**)) = *h*. This spectrum can be derived from the scaling behavior of partition functions defined in terms of WT coefficients (Muzy et al., 1991, 1994; Arneodo et al., 1995, 2008, 2003):

where *q* ∈ ℝ. Then, from the scaling function τ(*q*), *D*(*h*) is obtained by a Legendre transform:

In practice, to avoid instabilities in the estimation of the singularity spectrum *D*(*h*) through the Legendre transform, we instead compute the following expectation values (Muzy et al., 1994; Arneodo et al., 1995, 2003)

and

where $\hat{W}(q,l,a)=\frac{{M}{(\xb7,a)}^{q}}{Z(q,a)}$ corresponds to the Bolzmann weight in the analogy that connects the multifractal formalism to thermodynamics (Arneodo et al., 1995). Then, from the slopes of *h*(*q*, *a*) and *D*(*q*, *a*) vs. ln *a*, we get *h*(*q*) and *D*(*q*), and in turn the *D*(*h*) singularity spectrum as a curve parametrized by *q*.

### 2.4. Monofractal vs. Multifractal Functions

Homogeneous *monofractal* functions are functions with singularities of unique Hölder exponent *H*. Their τ(*q*) spectrum is a linear curve of slope *H*. Monofractal scaling indeed means that the shape of the probability distribution function (pdf) of rescaled wavelet coefficients (${M}$(·, *a*)/*a*^{H}) is independent of the scale *a*: ${\rho}_{{{M}}_{a}/{a}^{H}}(w)=\rho (w)$, where ρ(*w*) is a “universal” pdf (Arneodo et al., 2002, 2011). A nonlinear τ(*q*) is the signature of nonhomogeneous *multifractal* functions with Hölder exponent *h*(*t*) (resp. *h*(**x**)) fluctuating over time *t* (Muzy et al., 1991, 1994; Arneodo et al., 1995, 2008) [resp. over space **x** (Arneodo et al., 2000, 2003; Decoster et al., 2000; Roux et al., 2000)]. τ(*q*) data are generally well approximated by the log-normal quadratic spectrum:

where the coefficients *c*_{n} > 0. The corresponding singularity spectrum has a quadratic single-humped shape

where *c*_{0} = −τ(0) is the fractal dimension of the support of singularities of Σ (resp. *I*), *c*_{1} is the value of *h* that maximizes *D*(*h*) and the *intermittency coefficient c*_{2} (Delour et al., 2001) characterizes the width of the *D*(*h*) spectrum as an indication of a change in WT coefficient statistics across scales (Muzy et al., 1991, 1994; Arneodo et al., 1995, 2002, 2003, 2008).

## 3. Description of Data

### 3.1. Study Design and Population

Subjects were recruited for the present study from Perm Region Oncological Dispensary using procedures approved by the Local Ethics Committe (Gileva et al., 2012). Patients gave informed consent to participate in this study via recording IR thermograms and X-ray mammograms of their two mammary glands: the cancerous breast and the contralateral unaffected breast. The thermography (resp. mammography) database included 33 (resp. 30 among these 33) females, aged 37–83 (average 57 years), who all went through surgery to remove a histologically-confirmed malignant tumor (invasive ductal and/or lobular cancer) a few weeks after thermo-(mammo-)grams were recorded. The tumors were found at depths varying from 1 to 12 cm into the tissue, with a size (diameter) varying from 1.2 cm up to 6.5 cm (Table S1).

### 3.2. IR Thermography Imaging

Recording protocol was described in our preliminary study (Gerasimova et al., 2014). The patient's two breasts were imaged with an InSb photovoltaic (PV) detector (Joro et al., 2008). We performed imaging with the patient sitted with arms down to avoid too much discomfort. Frontal images were taken at a distance ~1 m from the patient, with an environmental room temperature ~20−22 °C. The image frame rate was 50 Hz. During the 10 min immobile imaging period, we collected 30,000 256 × 320 pixel^{2} image frames at 14 bits on a computer connected to the PV camera. We used skin surface markers as reference points for low-frequency motion correction in the analysis.

### 3.3. X-Ray Mammography

The mammographic procedure involved taking two X-ray images of the (two) compressed breasts of each patient, the standard medio-lateral oblique (MLO) and cranio-caudal (CC) views. Spatial resolution of the images was 50 μm per pixel. The ionizing radiation dose was 0.4 mSv per patient.

### 3.4. Data Sampling

#### 3.4.1. Thermograms

As commonly done for noise signals (Muzy et al., 1994; Audit et al., 2002), and previously performed to analyze rainfall time-series (Venugopal et al., 2006), the 1D WTMM method was applied to the cumulative (or integral) Σ of the temperature time-series (Figure S1), using the second-order compactly supported version ${\psi}_{(3)}^{(2)}$ of the Mexican hat (Roux et al., 1999) (see Figure S1 in (Gerasimova et al., 2014)). The singularities with possible negative Hölder exponent −1 < *h* < 0, became singularities with 0 < *h*_{c} = *h* + 1 < 1 in the cumulative. We grouped single-pixel temperature time-series (Figures S1A,B) into 8 × 8 pixel^{2} squares spanning 10 × 10 mm^{2} and covering the entire breast (see Figures 1A,A′). The results thus correspond to averaged partition functions and multifractal τ(*q*) (Equation 10) and *D*(*h*) (Equation 11) spectra computed from the WT skeleton (e.g., Figures S1E,F) of 64 cumulative temperature time-series in these 8 × 8 subareas. For each breast analyzed, we have color-coded these squares according to the diagnosis obtained with the 1D WTMM method (Table 2 in Gerasimova et al., 2014): monofractal (*c*_{2} < 0.03, red), multifractal (*c*_{2} ≥ 0.03, blue) and no scaling (white) (Figures 1A,A′).

**Figure 1. Wavelet-based multifractal segmentation of dynamic infrared thermograms and X-ray mammograms**. Patient 20 (age 56): cancerous right breast **(A–C)** and contralateral unaffected left breast **(A′–C′)**. **(A,A′)** As estimated from the τ(*q*) spectrum of skin temperature temporal fluctuations computed with the 1D WTMM method (Figures 2D–F), 8 × 8 pixel^{2} squares spanning 10 × 10 mm^{2} were color coded according to monofractal (*c*_{2} < 0.03, red), multifractal (*c*_{2} ≥ 0.03, blue) and no scaling (white) diagnostic, where *c*_{2} is the intermittency coefficient that defines the width of the *D*(*h*) singularity spectrum (Equations 10, 11) (Gerasimova et al., 2014). **(B,B′)** As estimated from the τ(*q*) spectra of CC mammographic view computed with the 2D WTMM method (Figures 2A–C), 256 × 256 pixel^{2} squares spanning 12.8 × 12.8 mm^{2} were color coded according to monofractal *H* < 0.45 (blue), 0.45 ≤ *H* ≤ 0.55 (yellow), *H* > 0.55 (red) and no scaling (pink). **(C,C′)** Same as **(B,B′)** for MLO mammographic view.

#### 3.4.2. Mammograms

For 2D WTMM analysis of mammograms, we used the isotropic Gaussian function

as a smoothing function, meaning that we worked with first-order (one vanishing moment) analyzing wavelets ψ_{1} and ψ_{2} (Equation 2). We divided the entire mammographic images into 360 × 360 pixels^{2} overlapping (to control edge effects) squares spanning 18 × 18 mm^{2}. Subimages of the breast in each of these squares were analyzed using the 2D WTMM methodology (Arneodo et al., 2000, 2003; Decoster et al., 2000; Roux et al., 2000). For a given square, when computing the partition functions and multifractal spectra from the WT skeleton, only the central (256 × 256 pixels^{2}) part that does not overlap with the central part of the neighboring squares was taken into account (Figure S2). For each breast analyzed, we have color-coded these squares according to the monofractal diagnosis obtained with the 2D WTMM method: long-range correlations (*H* > 0.55, red), anti-correlations (*H* < 0.45, blue), no-correlations (0.45 ≤ *H* ≤ 0.55, yellow), and no scaling (pink) (Figures 1B,B′,C,C′).

### 3.5. Statistical Tests

Statistical analyses were performed using the R statistical package (https://cran.r-project.org/). A nonparametric Wilcoxon rank-sum test was used to calculate *p*-values when comparing CC and MLO mammograms, cancer breasts and contralateral unaffected breasts, and thermograms and mammograms. Linear correlation of data sets were analyzed via a Pearson test.

### 3.6. Software and Documentation

The numerical procedure to perform the WTMM analysis of 1D signals is avaliable at http://perso.ens-lyon.fr/benjamin.audit/LastWave

LastWave is an open source software written in C. We recommend interested users read the LastWave C-Application Programming Interface documentation and to contact the corresponding author to be directed to the part of the code of most relevance to them. For 2D WTMM analysis, we adapted home-made Xsmurf software written in C, for specific application to mammograms.

## 4. Results

### 4.1. Mammographic Tissue Classification Using the 2D WTMM Method

We analyzed individual 360 × 360 pixel^{2} mammographic subimages of patient breasts with malignant tumor using the 2D WTMM method (Materials and Methods). From the central 256 × 256 pixels^{2} part of the WT skeleton we computed the partition functions *Z*(*q*, *a*) (Equation 6) that were found to display nice scaling properties over a range of scales 1 ≲ log_{2} *a* ≲ 3 corresponding to [0.7 mm, 2.8 mm] for linear regression fit estimates in a logarithmic representation (Figure 2A). The scaling deteriorates when considering larger scales due to finite size effects. In the range −1 ≲ *q* ≲ 3, statistical convergence is achieved; the so-obtained τ(*q*)-spectrum is remarkably linear, as the signature of monofractal roughness fluctuations and this wherever the location of the mammographic subimages over the cancerous breast (Figure 2B). Indeed the data are well-fitted by the theoretical spectrum τ(*q*) = *qH* − 2 of monofractal rough surfaces, which are almost everywhere singular with a unique hölder exponent *h* = *H* (Arneodo et al., 2000, 2003, 2008; Decoster et al., 2000). This is confirmed when computing the *D*(*h*) singularity spectrum from the slopes *h*(*q*) and *D*(*q*) of the partition functions *h*(*q*, *a*) (Equation 8) and *D*(*q*, *a*) (Equation 9) vs. ln *a* (Figure S3). *h*(*q*) and *D*(*q*) do not significantly depend on *q*, meaning that the *D*(*h*) singularity spectrum reduces to a single point *D*(*h* = *H*) = 2 (Figure 2C). Actually, what possibly changes when investigating different regions of cancerous breasts is the Hurst exponent *H* determined by the WTMM method (Figures 2B,C). As previously noted in a preliminary study (Kestener et al., 2001; Arneodo et al., 2003), we recovered *H* ≃ 1/3 in regions of fatty breast tissues as an indication of antipersistent (anti-correlated) roughness fluctuations, while dense tissues more likely displayed monofractal scaling with *H* ≃ 2/3 as an indication of persistent long-range correlations. But what was undetermined in previous analysis, was the existence of breast regions where the estimated Hurst exponent *H* ≃ 1/2, the hallmark of uncorrelated monofractal rough surfaces (Figures 2B,C). Importantly, the loss of correlations in roughness fluctuations is mainly observed in the region of the breast where the malignant tumor is located.

**Figure 2. Multifractal analysis of IR temperature time series and X-ray mammograms of patient 20**. Comparative analysis of spatial roughness fluctuations in individual 256 × 256 pixel^{2} squares in the cancerous right breast with respectively monofractal *H* < 0.45 (blue), 0.45 ≤ *H* ≤ 0.55 (yellow) and *H* > 0.55 (red): **(A)** log_{2} *Z*(*q*, *a*) vs. log_{2} *a*; **(B)** τ(*q*) vs. *q* estimated by linear regression fit of log_{2} *Z*(*q*, *a*) vs. log_{2} *a* (Equation 6) over a range of space-scales [0.7; 2.8] mm; **(C)** *D*(*h*) vs. *h*. Comparative analysis of temperature fluctuations in a 8 × 8 pixel^{2} square in the cancerous right breast (monofractal *c*_{2} < 0.03, red) and contralateral unaffected left breast (multifractal *c*_{2} ≥ 0.03, blue): **(D–F)** as **(A–C)**.

### 4.2. Segmentation of Breast Mammograms into Physiologically Altered (Risky) and Normal Regions: Comparative Analysis of the CC and MLO Mammographic Views

We extended our wavelet-based multifractal analysis of CC and MLO mammographic images to the entire sets of 256 × 256 pixels^{2} squares that cover every breast of the 30 (among 33) patients that proceeded through X-ray mammography prior to surgery. As shown in Figures 1B,B′,C,C′ for patient 20, we color-coded each of these squares according to the parameter *H* estimated with the 2D WTMM method: (i) *H* < 0.45, blue (“fatty”) anti-correlated rough surface, (ii) *H* > 0.55, red (“dense”) long-range correlated rough surface, (iii) 0.45 ≤ *H* ≤ 0.55, yellow uncorrelated rough surface, and pink when no convincing scaling was observed in the considered square. For the two mammographic views of each breast, we calculated the numbers of blue (*N*_{b}), red (*N*_{r}), yellow (*N*_{y}) and pink (*N*_{n}) squares and corresponding percentages of breast coverage (Table S2). We first confirmed that, except in a minority of squares, monofractal scaling was clear. For the cancerous breasts, the mean percentage of no-scaling pink squares was 6.5% (resp. 5.7%) in CC (resp. MLO) mammograms, which was consistent with the mean percentage 5.9% (resp. 6.1%) in CC (resp. MLO) mammograms obtained for the (a priori healthy) contralateral unaffected breasts. Most of the patients had breasts consisting of mostly fatty tissues, with a majority of anti-correlated blue squares covering on average 70.9% (resp. 72.4%) of CC (resp. MLO) mammograms of the cancerous breasts and 79.1% (resp. 80.1%) in CC (resp. MLO) mammograms of the contralateral unaffected breasts. This excess was compensated by a rather small percentage of long-range correlated red squares 6.4% (resp. 6.5%) in CC (resp. MLO) mammograms of the cancerous breasts, as well as in the CC (resp. MLO) mammograms of the contralateral unaffected breasts [6.5% (resp. 6.5%)]. A similar agreement was observed between the percentages of uncorrelated yellow squares (further called *H* = 0.5 squares) in the CC and corresponding MLO mammograms of each breast of each patient (Figure 3A and Table S2). However, the number (Figures S4, S5A,A′) and percentage (Figures 3A, 4A,A′) of uncorrelated *H* = 0.5 squares in the cancerous breast mammograms were both larger than in the contralateral unaffected breast mammograms. In the CC (resp. MLO) view, *N*_{y}/*N*_{Tot} ranges in the interval [4.5%, 41.4%] (resp. [4.9%, 32%]) with a mean percentage 16.1% (resp. 15.4%) for the breasts with the malignant tumor, as compared to the interval [0.6%, 27%] (resp. [0%, 16.7%]) and mean percentage 8.5% (resp. 7.3%) for the opposite breasts. A statistical comparison between cancerous breasts and contralateral unaffected breasts was performed using a Wilcoxon rank-sum test which yielded *p* < 10^{−4} and 10^{−5} for the moduli of the values in Figure 3A and Figure S4A, respectively. Moreover, the moduli of values for both cancerous breasts and contralateral unaffected breasts in CC and MLO mammograms exhibited linear correlation. In Figures 4A,A′, cancerous breasts and contralateral unaffected breasts were different with *p* < 10^{−8}. This statistically significant excess of *H* = 0.5 squares in the mammograms of cancerous breasts is evidence that the loss of correlations in the mammographic spatial density fluctuations can be used as an indicator of malignancy.

**Figure 3. Analysis of the percentage of monofractal uncorrelated H = 0.5 256 × 256 pixel^{2} squares in the mammograms of the breasts of 30 patients with breast cancer. (A)** Percentage of

*H*= 0.5 squares in cancerous (red) and contralateral unaffected (blue) breasts: MLO view vs. CC view.

**(B)**Mean percentage of

*H*= 0.5 squares in MLO and CC views: contralateral unaffected breast (CUB) vs. cancerous breast (CB). In

**(A,B)**, the numbers correspond to patient numbers defined in Table S1.

**Figure 4. Differential monofractal H = 0.5 signature on the (CC and MLO) mammograms of cancerous breasts (CB) and contralateral unaffected breasts (CUB). (A)** Histograms of the percentage of squares in the mammograms of the CBs (red) and CUBs (blue) of 30 patients with breast cancer.

**(B)**Histograms of the size of

*H*= 0.5 square clusters (see text) in CBs (red) and CUBs (blue).

**(C)**Histograms of the size of the largest

*H*= 0.5 square cluster in CBs (red) and CUBs (blue).

**(A′–C′)**same as

**(A–C)**for the corresponding cumulative distribution functions. Clusters are defined by squares sharing a common edge or corner. In

**(A–C)**, the pink and green vertical lines correspond to the mean (dashed line) and median (dashed-dotted line) of the histogram obtained for CBs and CUBs respectively.

### 4.3. Comparative Statistical Analysis of Mammogram Roughness Fluctuations in Cancerous and Contralateral Unaffected Breasts

Asymmetry of patient's breast can be a sign of breast cancer. When comparing the percentages of uncorrelated *H* = 0.5 squares on both the cancerous and contralateral unaffected breast mammograms of each patient, we found that a majority (25/30) have more *H* = 0.5 squares (average over CC and MLO views) on the cancerous breast (Figure 3B and Figure S4B and Table S2), including patients 20 (Figure 1), 17 (Figure S6), 24 (Figure S7), and 29 (Figure S8). For the 5 patients left including patients 3 (Figure S9), 7 (Figure 5) and 30 (Figure S10), an important and slightly larger (or equal) percentage of *H* = 0.5 squares was found in the second breast as a probable indication of some physiological and architectural changes in the undiagnosed opposite breast. When looking at the excess of *H* = 0.5 squares in the cancerous breasts we realized that these uncorrelated squares were not sparsely distributed all over the breast, but were concentrated in specific regions likely surrounding the underlying malignant tumor (Figures 1B,C and Figures S6B,C, S7B,C, S8B,C). To quantify this inhomogeneous distribution, we investigated the size distribution of *H* = 0.5 clusters defined as *H* = 0.5 squares sharing a common edge or corner (Table 1). If cancerous and contralateral unaffected breasts display about the same number of isolated *H* = 0.5 squares (singlet), the former are comprised of clusters with sizes larger than 2 and up to 26 (patient 13). The corresponding histogram of *H* = 0.5 cluster size displays a long tail in cancerous breasts that has no counterpart in the contralateral unaffected breasts (Figure 4B). As seen on the cumulative histogram (Figure 4B′), more than 20% of *H* = 0.5 clusters have a size larger than 5 in the cancerous breast as compared to only 8.7% in the contralateral unaffected breasts. This clustering effect is even more pronounced on the histogram of the size of the largest *H* = 0.5 square cluster (Figure 4C), with 33.3% of the cancerous breasts with a largest *H* = 0.5 square cluster bigger than 10 as compared to only 3.4% for the contralateral unaffected breasts (Figure 4C′). Using Wilcoxon rank-sum test to compare cancerous breasts and contralateral unaffected breasts in Figures 4A–C yielded *p* < 10^{−8}, 10^{−3}, and 10^{−8}, respectively. Note that the only contralateral unaffected breast with a *H* = 0.5 square cluster of size 22 is the second breast of patient 3 that was previously shown to have more *H* = 0.5 squares that the breast with the malignant tumor (Figure S9). Quite similar statistical results were obtained when defining *H* = 0.5 clusters from squares sharing a common edge only (Figures S5B,B′,C,C′ and Table S3).

**Figure 5. Wavelet-based multifractal segmentation of dynamic infrared thermograms and X-ray mammograms**. Same as Figure 1 but for patient 7 (age 47): cancerous right breast **(A–C)** and contralateral unaffected left breast **(A′–C′)**.

**Table 1. Statistics of monofractal uncorrelated H = 0.5 yellow square clusters in the CC and MLO mammographic views of the cancerous and contralateral unaffected breasts of our patients with breast cancer**.

### 4.4. Comparative Multifractal Analysis of the Skin Temperature Dynamics and Mammogram Spatial Roughness Fluctuations

For the same cohort of patients, the analysis of dynamic IR thermograms put into light some significant change in the nature of temporal fluctuations of breast skin temperature in cancerous breasts (Gerasimova et al., 2013, 2014). When using the 1D WTMM method, the computation of the partition functions *Z*(*q*, *a*) (Equation 6), *h*(*q*, *a*) (Equation 8) and *D*(*q*, *a*) (Equation 9) revealed that, besides the cardiogenic and vasomotor perfusion oscillations, the temperature dynamics display some scale invariant background component that extends from the characteristic human respiratory frequency (≳ 0.3 Hz) up to the cross-over frequency (≲ 4 Hz) toward instrumental white noise (Figure 2D). For normal breasts as well as in healthy regions of cancerous breasts, the τ(*q*) spectrum is definitely non-linear (Figure 2E) and accordingly the *D*(*h*) spectrum has a single humped shape (Figure 2F), the hallmark of multifractal scaling. By contrast, in the region of the malignant tumor in cancerous breasts the τ(*q*) spectrum is nearly linear and the *D*(*h*) spectrum reduces to a single point as the signature of monofractal scaling. As illustrated in Figures 1A,A′ for patient 20, we have segmented the two breasts of each patient into 8 × 8 pixels squares that were further color-coded according to the monofractal (red), multifractal (blue) or no-scaling (white) diagnosis obtained with the 1D WTMM method. In good agreement with our previous segmentation of CC (Figures 1B,B′) and MLO (Figures 1C,C′) mammograms, more “risky” monofractal squares (49.7% coverage) were found in the cancerous breast (Figure 1A) than in the contralateral unaffected breast (7.7%) (Figure 1A′). In addition to their number, the location and clustering of these “risky” squares, where the multifractal complexity of temperature fluctuations is lost, correlate with the clustering of the “risky” *H* = 0.5 squares in CC and MLO mammograms indicating some additional evidence of loss of correlations in spatial density fluctuations in the upper outer quadrant of the right breast of patient 20 where the underlying tumor is located.

The results obtained for patient 20 are quite representative of the outcome of our overall comparative analysis of IR thermograms and X-ray mammograms for our set of 30 patients (Figure 6, Figure S11). The use of Wilcoxon rank-sum test to compare cancerous breasts and contralateral unaffected breasts in Figure 6 and Figure S11 yielded *p* < 10^{−4} and 10^{−3}, respectively. But both types of breasts did not exhibit linear correlation between mammogram and thermogram data. Most patients have consistent large numbers (breast coverage ≳ 10%) of “risky” uncorrelated *H* = 0.5 (yellow) squares in mammograms and monofractal (red) squares in thermograms of their cancerous breast as compared to their contralateral unaffected breast (breast coverage ≲ 10%), including patients 17 (Figure S6) and 29 (Figure S8). Note that for the subset of (5) patients, including patients 3 (Figure S9), 7 (Figure 5), and 30 (Figure S10), who were shown to have also a lot of “risky” *H* = 0.5 squares on their contralateral unaffected breast mammograms, they also have a lot of risky monofractal squares on the thermograms of their two breasts which might be an additional sign of the possible extension of cancer to their second breast. Among the patients with few “risky” monofractal squares (≲ 10%) on their IR thermograms of their cancerous breast, 4 correspond to rather deep tumors, namely patients 12 (size 1.8 cm, depth 12 cm), 16 (3.4 cm, 7 cm), 18 (3.49 cm, 6 cm) (Figure 7) and 28 (3.49 cm, 8 cm). For these 4 patients the percentage of “risky” *H* = 0.5 squares on their corresponding cancerous breast mammograms is significantly high (≳ 10%), meaning that, although imperceptible in temperature dynamics at the skin surface, the change in the microenvironment of these deep tumors turns out to be detectable with X-ray mammography. Let us also mention that if the Paget disease of the nipple of the left breast of patient 24 (Figure S7) was hardly identified in the IR thermograms with only few (breast coverage 4.9%), but well localized, “risky” monofractal squares (Figure S7A), it could not be missed on the CC (Figure S7B) and MLO (Figure S7C) mammograms, with respective 41.4% and 24% coverages by “risky” uncorrelated *H* = 0.5 squares.

**Figure 6. Comparative analysis of thermograms and (CC and MLO averaged) mammograms**. Percentage of monofractal squares in thermograms vs. percentage of *H* = 0.5 squares in mammograms of the cancerous (red) and contralateral unaffected (blue) breasts of the 30 patients with breast cancer. The numbers correspond to patient numbers defined in Table S1.

**Figure 7. Wavelet-based multifractal segmentation of dynamic infrared thermograms and X-ray mammograms**. Same as Figure 1 but for patient 18 (age 37): cancerous left breast **(A–C)** and contralateral unaffected right breast **(A′–C′)**.

### 4.5. Comparative Analysis of the Mammograms of a Cancerous Breast before and 11 Months after Surgery: a Control Experiment

As a control, we reproduced our 2D WTMM multifractal analysis on the mammograms of the two breasts of patient 30 that were recorded 11 months after surgical treatment of invasive ductal carcinoma of the right breast. For both mammographic views, the number of “risky” *H* = 0.5 squares in the cancerous breast was reduced from *N*_{y} = 8 (Figure 8A) to 3 (Figure 8C) in CC mammograms, and from *N*_{y} = 13 (Figure 8B) to 1 (Figure 8D) in MLO mammograms. However, less expected was the quite similar quantitative reduction observed in the other breast of patient 30, who was one of the 5 patients having slightly more *H* = 0.5 squares on the contralateral unaffected breast than on the cancerous breast (Figure S10): *N*_{y} = 12 (Figure 8A′) to 0 (Figure 8C′) in CC mammograms, and *N*_{y} = 9 (Figure 8B′) to 3 (Figure 8D′) in MLO mammograms. These observations suggest that the surgery and associated therapeutic treatment were efficient with regards to removing the malignant tumor, as well as of its tumorigenic microenvironment including cancer stem cells. They also attest to some regeneration of the breast tissue, probably via the activation of quiescent epithelial stem cells housed in the ducts, and which were shown to exhibit a high proliferative, self-renewal and morphogenic capacity in culture (Villadsen et al., 2007). The underlying recovery mechanisms appear to be engaged and operating not only in the cancerous breast but also in the other breast which contained an above average number of “risky” *H* = 0.5 squares on both mammographic views.

**Figure 8. Comparative multifractal analysis of the mammograms of the cancerous breast of patient 30 (age 54) before surgery (A,B) and 11 months after surgery (C,D)**. The color coding of the 256 × 256 pixel^{2} squares in the CC **(A,C)** and MLO **(B,D)** mammogram views is the same as in Figures 1B,B′,C,C′. **(A′–D′)** Similar comparative multifractal analysis of the mammograms of the contralateral unaffected breast of patient 30.

## 5. Conclusions

To summarize, we showed that the wavelet-based multifractal analysis of X-ray mammograms was able to detect the presence of an internal malignant tumor via some loss of correlations in the breast density spatial fluctuations. We demonstrated that this important change in mammogram roughness fluctuations strongly correlated with some drastic simplification of temperature dynamics recorded with an IR camera, from multifractal to homogeneous monofractal temporal fluctuations. The nature of this study was exploratory, with a data set limited to females who all went through surgery to remove the histologically confirmed malignant tumor of rather large size (diameter between 1.2 and 6.5 cm). To our knowledge, our study is the first to report on the observation of some physiological alteration and architectural disorganization of the microenvironmet of breast tumors using classical screening techniques, including the currently, most used X-ray mammography. Since the alteration of the niche surrounding a breast tumor is likely to correspond to a transition from an antitumorigenic to a tumorigenic microenvironment that probably preceeds and further facilitates the process of oncogenic transformation and tumor progression (Bissell et al., 2005; Rønnov-Jessen and Bissell, 2009; Bissell and Hines, 2011; Maguer-Satta, 2011; Lu et al., 2012), the results reported are not only novel, but show great promise toward improving early breast cancer diagnosis that is known to be critical for the treatment and survival of the patient. Of course these results deserve to be confirmed over a larger set of patients at different stages of cancer development. We just started the analysis of longitudinal studies that should bring statistical estimates of the sensitivity and specificity of our mathematical and computational approach mainly based on the estimate of the mammographic density fluctuation index *H*. This preliminary study does suggest that combining sparse and sometimes painful and quite uncomfortable mammography examinations with more frequent inexpensive, quick and painless IR thermography examinations could become an efficient routine breast cancer screening method to identify, as early as possible, women with high risk of breast cancer.

## Author Contributions

Conception and design: AA, ON, AK, FA, OG. Development and methodology: AA, ON, EG, BA, AK. Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): EG. Analysis and interpretation of data: AA, EG, BT, ZM, BA, SR, AK, FA, OG, ON. Writing, review, and/or revision of the manuscript: AA. Administrative, technical or material support (i.e., requiling and organizing data, constructing databases): EG, BT, ZM, BA, FA. Study supervision: AA.

## Funding

This work was supported by INSERM, ITMO Cancer for its financial support under contract PC201201-084862 “Physiques, mathématiques ou sciences de l'ingénieur appliqués au Cancer,” the Russian Foundation for Basic Research (grant 16-41-590235), the Perm Regional Government (Russia) with the contract “Multiscale approaches in mechanobiology for early cancer diagnosis” and the Maine Cancer Foundation.

## Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

## Acknowledgments

We are thankful to Prof. O. Orlov for his support in IR thermal and X-ray mammographic imaging at Perm Region Cancer Hospital and to K. A. Batchelder, Y. Bayandin, I. Panteleev, O. Plekhov, and S. Uvarov for fruitful discussions.

## Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/article/10.3389/fphys.2016.00336

## References

Antoine, J. P., Murenzi, R., Vandergheynst, P., and Ali, S. T. (2008). *Two-Dimensional Wavelets and their Relatives*. Cambridge: Cambridge University Press.

Arneodo, A., Audit, B., Decoster, N., Muzy, J. F., and Vaillant, C. (2002). “A wavelet based multifractal formalism: application to DNA sequences, satellite images of the cloud structure and stock market data,” in *The Science of Disasters: Climate Disruptions, Heart Attacks, and Market Crashes*, eds A. Bunde, J. Kropp, and H. J. Schellnhuber (Berlin: Springer Verlag), 26–102.

Arneodo, A., Audit, B., Kestener, P., and Roux, S. G. (2008). Wavelet-based multifractal analysis. *Scholarpedia* 3:4103. doi: 10.4249/scholarpedia.4103

Arneodo, A., Bacry, E., and Muzy, J.-F. (1995). The thermodynamics of fractals revisited with wavelets. *Physica A* 213, 232–275. doi: 10.1016/0378-4371(94)00163-N

Arneodo, A., Decoster, N., Kestener, P., and Roux, S. G. (2003). A wavelet-based method for multifractal image analysis: from theoretical concepts to experimental applications. *Adv. Imaging Electr. Phys.* 126, 1–92. doi: 10.1016/S1076-5670(03)80014-9

Arneodo, A., Decoster, N., and Roux, S. G. (2000). A wavelet-based method for multifractal image analysis. I. Methodology and test applications on isotropic and anisotropic random rough surfaces. *Eur. Phys. J.* 15, 567–600. doi: 10.1007/s100510051161

Arneodo, A., Grasseau, G., and Holschneider, M. (1988). Wavelet transform of multifractals. *Phys. Rev. Lett.* 61, 2281–2284. doi: 10.1103/PhysRevLett.61.2281

Arneodo, A., Grasseau, G., and Kostelich, E. J. (1987). Fractal dimensions and *f*(α) spectrum of the Hénon attractor. *Phys. Lett. A* 124, 426–432. doi: 10.1016/0375-9601(87)90546-9

Arneodo, A., Vaillant, C., Audit, B., Argoul, F., d'Aubenton-Carafa, Y., and Thermes, C. (2011). Multi-scale coding of genomic information: from DNA sequence to genome structure and function. *Phys. Rep.* 498, 45–188. doi: 10.1016/j.physrep.2010.10.001

Audit, B., Bacry, E., Muzy, J.-F., and Arneodo, A. (2002). Wavelet-based estimators of scaling behavior. *IEEE Trans. Inform. Theory* 48, 2938–2954. doi: 10.1109/TIT.2002.802631

Audit, B., Baker, A., Chen, C.-L., Rappailles, A., Guilbaud, G., Julienne, H., et al. (2013). Multiscale analysis of genome-wide replication timing profiles using a wavelet-based signal-processing algorithm. *Nat. Protoc.* 8, 98–110. doi: 10.1038/nprot.2012.145

Ayer, T., Ayvaci, M. U., Liu, Z. X., Alagoz, O., and Burnside, E. S. (2010). Computer-aided diagnostic models in breast cancer screening. *Imaging Med.* 2, 313–323. doi: 10.2217/iim.10.24

Bacry, E., Muzy, J. F., and Arneodo, A. (1993). Singularity spectrum of fractal signals from wavelet analysis: exact results. *J. Stat. Phys.* 70, 635–674. doi: 10.1007/BF01053588

Batchelder, K. A., Tanenbaum, A. B., Albert, S., Guimond, L., Kestener, P., Arneodo, A., et al. (2014). Wavelet-based 3D reconstruction of microcalcification clusters from two mammographic views: new evidence that fractal tumors are malignant and euclidean tumors are benign. *PLoS ONE* 9:e107580. doi: 10.1371/journal.pone.0107580

Besançon, R., Valsesia-Wittmann, S., Puisieux, A., de Fromentel, C. C., and Maguer-Satta, V. (2009). Cancer stem cells: the emerging challenge of drug targeting. *Curr. Med. Chem.* 16, 394–416. doi: 10.2174/092986709787315531

Bissell, M. J., and Hines, W. C. (2011). Why don't we get more cancer? A proposed role of the microenvironment in restraining cancer progression. *Nat. Med.* 17, 320–329. doi: 10.1038/nm.2328

Bissell, M. J., Kenny, P. A., and Radisky, D. C. (2005). Microenvironmental regulators of tissue structure and function also regulate tumor induction and progression: the role of extracellular matrix and its degrading enzymes. *Cold Spring Harb. Symp. Quant. Biol.* 70, 343–356. doi: 10.1101/sqb.2005.70.013

Bissell, M. J., and Labarge, M. A. (2005). Context, tissue plasticity, and cancer: Are tumor stem cells also regulated by the microenvironment? *Cancer Cell* 7, 17–23. doi: 10.1016/j.ccr.2004.12.013

Decoster, N., Roux, S. G., and Arneodo, A. (2000). A wavelet-based method for multifractal image analysis. II. Applications to synthetic multifractal rough surfaces. *Eur. Phys. J.* 15, 739–764. doi: 10.1007/s100510051179

Delour, J., Muzy, J.-F., and Arneodo, A. (2001). Intermittency of 1D velocity spatial profiles in turbulence: a magnitude cumulant analysis. *Eur. Phys. J. B* 23, 243–248. doi: 10.1007/s100510170074

Demicheli, R. (2001). Tumour dormancy: findings and hypotheses from clinical research on breast cancer. *Semin. Cancer Biol.* 11, 297–306. doi: 10.1006/scbi.2001.0385

Faraldo, M. M., Teulière, J., Deugnier, M. A., Taddei-De La Hosseraye, I., Thiery, J. P., and Glukhova, M. A. (2005). Myoepithelial cells in the control of mammary development and tumorigenesis: data from genetically modified mice. *J. Mammary Gland Biol. Neoplasia* 10, 211–219. doi: 10.1007/s10911-005-9582-8

Fenton, J. J., Abraham, L., Taplin, S. H., Geller, B. M., Carney, P. A., D'Orsi, C., et al. (2011). Effectiveness of computer-aided detection in community mammography practice. *J. Natl. Cancer Inst.* 103, 1152–1161.

Flynn, C. M., and Kaufman, D. S. (2007). Donor cell leukemia: insight into cancer stem cells and the stem cell niche. *Blood* 109, 2688–2692. doi: 10.1182/blood-2006-07-021980

Fuchs, E., Tumbar, T., and Guasch, G. (2004). Socializing with the neighbors: stem cells and their niche. *Cell* 116, 769–778. doi: 10.1016/S0092-8674(04)00255-7

Ganesan, K., Acharya, U. R., Chua, C. K., Min, L. C., Abraham, K. T., and Ng, K. H. (2013). Computer-aided breast cancer detection using mammograms: a review. *IEEE Rev. Biomed. Eng.* 6, 77–98. doi: 10.1109/RBME.2012.2232289

Gerasimova, E., Audit, B., Roux, S., Khalil, A., Argoul, F., Naimark, O., et al. (2013). Multifractal analysis of dynamic infrared imaging of breast cancer. *Europhys. Lett.* 104:68001. doi: 10.1209/0295-5075/104/68001

Gerasimova, E., Audit, B., Roux, S. G., Khalil, A., Gileva, O., Argoul, F., et al. (2014). Wavelet-based multifractal analysis of dynamic infrared thermograms to assist in early breast cancer diagnosis. *Front. Physiol.* 5:176. doi: 10.3389/fphys.2014.00176

Gileva, O. S., Freynd, G. G., Orlov, O. A., Libik, T. V., Gerasimova, E. I., Plekhov, O. A., et al. (2012). Interdisciplinary approaches to early diagnosis and screening of tumors and precancerous diseases (for example breast cancer). *RFBR J.* 74–75, 93–99.

Goldberger, A. L., Amaral, L. A. N., Hausdorff, J. M., Ivanov, P. C., Peng, C.-K., and Stanley, H. E. (2002). Fractal dynamics in physiology: alterations with disease and aging. *Proc. Natl. Acad. Sci. U.S.A.* 99, 2466–2472. doi: 10.1073/pnas.012579499

Goody, M. F., Kelly, M. W., Lessard, K. N., Khalil, A., and Henry, C. A. (2010). Nrk2b-mediated NAD+ production regulates cell adhesion and is required for muscle morphogenesis *in vivo*: Nrk2b and NAD+ in muscle morphogenesis. *Dev. Biol.* 344, 809–826. doi: 10.1016/j.ydbio.2010.05.513

Grant, J., Verrill, C., Coustham, V., Arneodo, A., Palladino, F., Monier, K., et al. (2010). Perinuclear distribution of heterochromatin in developing C. elegans embryos. *Chrom. Res.* 18, 873–885. doi: 10.1007/s10577-010-9175-2

Häberle, L., Wagner, F., Fasching, P. A., Jud, S. M., Heusinger, K., Loehberg, C. R., et al. (2012). Characterizing mammographic images by using generic texture features. *Breast Cancer Res.* 14:R59. doi: 10.1186/bcr3163

Hall, B., Andreeff, M., and Marini, F. (2007). The participation of mesenchymal stem cells in tumor stroma formation and their application as targeted-gene delivery vehicles. *Handb. Exp. Pharmacol.* 180, 263–283. doi: 10.1007/978-3-540-68976-8_12

Ho, A. D., and Wagner, W. (2007). The beauty of asymmetry: asymmetric divisions and self-renewal in the haematopoietic system. *Curr. Opin. Hematol.* 14, 330–336. doi: 10.1097/MOH.0b013e3281900f12

Ivanov, P. C., Amaral, L. A. N., Goldberger, A. L., Havlin, S., Rosenblum, M. G., Stanley, H. E., et al. (2001). From 1/f noise to multifractal cascades in heartbeat dynamics. *Chaos* 11, 641–652. doi: 10.1063/1.1395631

Ivanov, P. C., Amaral, L. A. N., Goldberger, A. L., Havlin, S., Rosenblum, M. G., Struzik, Z. R., et al. (1999). Multifractality in human heartbeat dynamics. *Nature* 399, 461–465. doi: 10.1038/20924

Jalalian, A., Mashohor, S. B. T., Mahmud, H. R., Saripan, M. I. B., Ramli, A. R. B., and Karasfi, B. (2013). Computer-aided detection/diagnosis of breast cancer in mammography and ultrasound: a review. *Clin. Imaging* 37, 420–426. doi: 10.1016/j.clinimag.2012.09.024

Jørgensen, K. J., and Gøtzsche, P. C. (2009). Overdiagnosis in publicly organised mammography screening programmes: systematic review of incidence trends. *BMJ* 339:b2587. doi: 10.1136/bmj.b2587

Joro, R., Lääperi, A.-L., Dastidar, P., Soimakallio, S., Kuukasjäri, T., Toivonen, T., et al. (2008). Imaging of breast cancer with mid- and long-wave infrared camera. *J. Med. Eng. Technol.* 32, 189–197. doi: 10.1080/03091900701234358

Karahaliou, A. N., Boniatis, I. S., Skiadopoulos, S. G., Sakellaropoulos, F. N., Arikidis, N. S., Likaki, L. E. A., et al. (2008). Breast cancer diagnosis: analyzing texture of tissue surrounding microcalcifications. *IEEE Trans. Inf. Technol. Biomed.* 12, 731–738. doi: 10.1109/TITB.2008.920634

Kestener, P., and Arneodo, A. (2003). Three-dimensional wavelet-based multifractal method: the need for revisiting the multifractal description of turbulence dissipation data. *Phys. Rev. Lett.* 91:194501. doi: 10.1103/PhysRevLett.91.194501

Kestener, P., and Arneodo, A. (2004). Generalizing the wavelet-based multifractal formalism to random vector fields: application to three-dimensional turbulence velocity and vorticity data. *Phys. Rev. Lett.* 93:044501. doi: 10.1103/PhysRevLett.93.044501

Kestener, P., Lina, J.-M., Saint-Jean, P., and Arneodo, A. (2001). Wavelet-based multifractal formalism to assist in diagnosis in digitized mammograms. *Image Anal. Stereol.* 20, 169–174. doi: 10.5566/ias.v20.p169-174

Khalil, A., Aponte, C., Zhang, R., Davisson, T., Dickey, I., Engelman, D., et al. (2009). Image analysis of soft-tissue in-growth and attachment into highly porous alumina ceramic foam metals. *Med. Eng. Phys.* 31, 775–783. doi: 10.1016/j.medengphy.2009.02.007

Khalil, A., Grant, J. L., Caddle, L. B., Atzema, E., Mills, K. D., and Arneodo, A. (2007). Chromosome territories have a highly nonspherical morphology and nonrandom positioning. *Chrom. Res.* 15, 899–916. doi: 10.1007/s10577-007-1172-8

Lee, C. H. (2002). Screening mammography: proven benefit, continued controversy. *Radiol. Clin. N. Am.* 40, 395–407. doi: 10.1016/S0033-8389(01)00015-X

Lee, J. T., and Herlyn, M. (2007). Microenvironmental influences in melanoma progression. *J. Cell. Biochem.* 101, 862–872. doi: 10.1002/jcb.21204

Li, L., and Neaves, W. B. (2006). Normal stem cells and cancer stem cells: the niche matters. *Cancer Res.* 66, 4553–4557. doi: 10.1158/0008-5472.CAN-05-3986

Li, Z., and Li, L. (2006). Understanding hematopoietic stem-cell microenvironments. *Trends Biochem. Sci.* 31, 589–595. doi: 10.1016/j.tibs.2006.08.001

Lu, P., Weaver, V. M., and Werb, Z. (2012). The extracellular matrix: a dynamic niche in cancer progression. *J. Cell Biol.* 196, 395–406. doi: 10.1083/jcb.201102147

Maguer-Satta, V. (2011). “The stem cell niche: the black master of cancer,” in *Cancer Stem Cells Theories and Practice*, ed S. Shostak (Rijeka: InTech). doi: 10.5772/13622

Marthiens, V., Kazanis, I., Moss, L., Long, K., and Ffrench-Constant, C. (2010). Adhesion molecules in the stem cell niche – more than just staying in shape? *J. Cell Sci.* 123, 1613–1622. doi: 10.1242/jcs.054312

Meselhy Eltoukhy, M., Faye, I., and Belhaouari Samir, B. (2012). A statistical based feature extraction method for breast cancer diagnosis in digital mammogram using multiresolution representation. *Comput. Biol. Med.* 42, 123–128. doi: 10.1242/jcs.054312

Moore, K. A., and Lemischka, I. R. (2006). Stem cells and their niches. *Science* 311, 1880–1885. doi: 10.1126/science.1110542

Morrison, S. J., and Kimble, J. (2006). Asymmetric and symmetric stem-cell divisions in development and cancer. *Nature* 441, 1068–1074. doi: 10.1038/nature04956

Morrison, S. J., and Spradling, A. C. (2008). Stem cells and niches: mechanisms that promote stem cell maintenance throughout life. *Cell* 132, 598–611. doi: 10.1016/j.cell.2008.01.038

Muzy, J.-F., Bacry, E., and Arneodo, A. (1991). Wavelets and multifractal formalism for singular signals: application to turbulence data. *Phys. Rev. Lett.* 67, 3515–3518. doi: 10.1103/PhysRevLett.67.3515

Muzy, J.-F., Bacry, E., and Arneodo, A. (1993). Multifractal formalism for fractal signals: the structure-function approach versus the wavelet-transform modulus-maxima methods. *Phys. Rev. E* 47, 875–884. doi: 10.1103/PhysRevE.47.875

Muzy, J.-F., Bacry, E., and Arneodo, A. (1994). The multifractal formalism revisited with wavelets. *Int. J. Bifurc. Chaos* 4, 245–302. doi: 10.1142/S0218127494000204

Nass, S. J., Henderson, C. I., and Lashof, J. C. (2001). *Mammography and Beyond : Developing Technologies for the Early Detection of Breast Cancer*. Washington, DC: National Academy Press.

Nicolay, S., Brodie of Brodie, E. B., Touchon, M., Audit, B., d'Aubenton-Carafa, Y., Thermes, C., et al. (2007). Bifractality of human DNA strand-asymmetry profiles results from transcription. *Phys. Rev. E* 75:032902. doi: 10.1103/PhysRevE.75.032902

Richard, C. D., Tanenbaum, A., Audit, B., Arneodo, A., Khalil, A., and Frankel, W. N. (2015). Swdreader: a wavelet-based algorithm using spectral phase to characterize spike-wave morphological variation in genetic models of absence epilepsy. *J. Neurosci. Methods* 242, 127–140. doi: 10.1016/j.jneumeth.2014.12.016

Rønnov-Jessen, L., and Bissell, M. J. (2009). Breast cancer by proxy: can the microenvironment be both the cause and consequence? *Trends Mol. Med.* 15, 5–13. doi: 10.1016/j.molmed.2008.11.001

Roux, S., Muzy, J.-F., and Arneodo, A. (1999). Detecting vorticity filaments using wavelet analysis: about the statistical contribution of vorticity filaments to intermittency in swirling turbulent flows. *Eur. Phys. J. B* 8, 301–322. doi: 10.1007/s100510050694

Roux, S. G., Arneodo, A., and Decoster, N. (2000). A wavelet-based method for multifractal image analysis. III. Applications to high-resolution satellite images of cloud structure. *Eur. Phys. J.* 15, 765–786. doi: 10.1007/s100510051180

Siegel, R. L., Miller, K. D., and Jemal, A. (2015). Cancer statistics, 2015. *CA Cancer J. Clin.* 65, 5–29. doi: 10.3322/caac.21254

Snow, C. J., Goody, M., Kelly, M. W., Oster, E. C., Jones, R., Khalil, A., et al. (2008). Time-lapse analysis and mathematical characterization elucidate novel mechanisms underlying muscle morphogenesis. *PLoS Genet.* 4:e1000219. doi: 10.1371/journal.pgen.1000219

Trumpp, A., Essers, M., and Wilson, A. (2010). Awakening dormant haematopoietic stem cells. *Nat. Rev. Immunol.* 10, 201–209. doi: 10.1038/nri2726

Tsai, N.-C., Chen, H.-W., and Hsu, S.-L. (2011). Computer-aided diagnosis for early-stage breast cancer by using Wavelet Transform. *Comp. Med. Imaging Graph.* 35, 1–8. doi: 10.1016/j.compmedimag.2010.08.005

Tysnes, B. B., and Bjerkvig, R. (2007). Cancer initiation and progression: involvement of stem cells and the microenvironment. *Biochim. Biophys. Acta* 1775, 283–297. doi: 10.1016/j.bbcan.2007.01.001

Venugopal, V., Roux, S. G., Foufoula-Georgiou, E., and Arneodo, A. (2006). Revisiting multifractality of high-resolution temporal rainfall using a wavelet-based formalism. *Water Resour. Res.* 42:W06D14. doi: 10.1029/2005WR004489

Keywords: breast cancer, dynamic IR thermograms, X-ray mammograms, multifractal analysis, wavelet transform modulus maxima method

Citation: Gerasimova-Chechkina E, Toner B, Marin Z, Audit B, Roux SG, Argoul F, Khalil A, Gileva O, Naimark O and Arneodo A (2016) Comparative Multifractal Analysis of Dynamic Infrared Thermograms and X-Ray Mammograms Enlightens Changes in the Environment of Malignant Tumors. *Front. Physiol*. 7:336. doi: 10.3389/fphys.2016.00336

Received: 06 June 2016; Accepted: 20 July 2016;

Published: 09 August 2016.

Edited by:

Plamen Ch. Ivanov, Harvard Medical School, USACopyright © 2016 Gerasimova-Chechkina, Toner, Marin, Audit, Roux, Argoul, Khalil, Gileva, Naimark and Arneodo. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) or licensor are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Alain Arneodo, alain.arneodo@u-bordeaux.fr