Shannon Entropy-Based Wavelet Transform Method for Autonomous Coherent Structure Identification in Fluid Flow Field Data

The coherent secondary flow structures (i.e., swirling motions) in a curved artery model possess a variety of spatio-temporal morphologies and can be encoded over an infinitely-wide range of wavelet scales. Wavelet analysis was applied to the following vorticity fields: (i) a numerically-generated system of Oseen-type vortices for which the theoretical solution is known, used for bench marking and evaluation of the technique; and (ii) experimental two-dimensional, particle image velocimetry data. The mother wavelet, a two-dimensional Ricker wavelet, can be dilated to infinitely large or infinitesimally small scales. We approached the problem of coherent structure detection by means of continuous wavelet transform (CWT) and decomposition (or Shannon) entropy. The main conclusion of this study is that the encoding of coherent secondary flow structures can be achieved by an optimal number of binary digits (or bits) corresponding to an optimal wavelet scale. The optimal wavelet-scale search was driven by a decomposition entropy-based algorithmic approach and led to a threshold-free coherent structure detection method. The method presented in this paper was successfully utilized in the detection of secondary flow structures in three clinically-relevant blood flow scenarios involving the curved artery model under a carotid artery-inspired, pulsatile inflow condition. These scenarios were: (i) a clean curved artery; (ii) stent-implanted curved artery; and (iii) an idealized Type IV stent fracture within the curved artery. Entropy 2015, 17 6618


Introduction
The definition of coherent (vortical) structure by Robinson [1] states that with reference to a two-dimensional velocity vector field, usually associated with regions of local swirling motion, in a frame of reference moving with the convection velocity, a vortex might be recognized as a region where streamlines describe spiral patterns.Such a characterization of the organized vortical structures, often called coherent structures, is commonly associated with regions of concentrated vorticity [2,3].
The identification of hemodynamic secondary flow (vortical) structures can be performed by the user probing the experimental data, applying concepts of flow physics and critical point theory [4][5][6][7][8].There are several methods of coherent structure detection that require the computation of the velocity gradient tensor [9].The Q-criterion defines a vortex as a spatial region where the Euclidean norm of the vorticity tensor dominates that of the rate of strain [4,10].The λ ci -criterion is a vortex identification method performed by critical-point analysis of the local velocity gradient tensor and its corresponding eigenvalues [5,7,10].All of these methods find patterns of swirling flow that include vortices with a high circulation and circular morphologies and strain-dominated flows with low circulation and elongated morphologies.However, the described methodologies provide a scale decomposition (in space or time) of the velocity vector field, all based on the analysis of the filtered counterpart of the raw dataset, wherein the frequency cut-off (or threshold) has to be tuned appropriately; the process involves subjectivity and manual intervention.Furthermore, for weakly turbulent flows, the typical wavenumber (cut-off frequency or threshold) of coherent structures is not known a priori, and hence, without any a priori assumptions, such thresholds cannot be established [3].
Longo [11] elucidated that Fourier transform-based flow analysis using harmonic sine-cosine functions space as bases produces the contribution of localized structures without carrying any information on their occurrence in the time (or space) series.Wavelet decomposition methods provide an advantage of resolving the total velocity field scale-by-scale, since coherent structures do not have a single preferential length scale, even in transitional flow regimes where vortices tend to cascade from large to smaller eddies [3,12,13].Especially in turbulent flows, the cascade models rely on that very assumption that energy density is segmented in the wavenumber domain and is exchanged between eddies of individual eddies, making wavelet decomposition methods have the potential to reveal such eddies [11].Exhaustive literature on the development of such methodologies aimed at the identification and characterization of coherent structures have been provided by Camussi and Guj [2], Camussi [3], Longo [11], Bonnet et al. [14] and Longo [15].
The main focus of this paper is to introduce an autonomous, threshold-free method of coherent structure identification.The keys ideas are built on the wavelet property of localization and extended using entropy-aided signal processing methodologies.The purpose of the study is to present an algorithmic approach using the tenets of the classical theory of communication, wherein experimental data or "signals" are acquired from non-invasive, hydrodynamic experiments that were geared toward investigation of hemodynamic secondary flows in a curved artery model.The ideas of entropy of the signal or signal uncertainty and its measure (bits per signal units) as defined in communication theory are retained [16][17][18].The outcome of the algorithmic approach presented herein is the realization of an optimal wavelet scale that resolved the coherent structures in two-dimensional planar regions of the curved artery model.
Vortical flow patterns that arise due to a combination of centrifugal instabilities and adverse pressure gradients in curved channels and tubes are termed secondary flow structures in the hydrodynamics literature.In clinical studies, they are referred as "spiral blood flow structures".Their prevalence in the cardiovasculature is as ubiquitous as the curvatures in the arterial network itself.They have the potential to affect the hemodynamics of the aortic arch, play a protective role toward arterial wall damage and may possess beneficial effects on the mechanisms of endothelial function.The clinical evidence of these structures was elucidated by Sengupta et al. [19], Stonebridge et al. [20] and Stonebridge and Brophy [21].In addition, these secondary flow structures are also common for a variety of applications ranging from heat and mass transfer enhancements to curved channel flows [13].
The belief that vortical flows and vortex-induced phenomena in the blood flows have a role to play in the proliferation of several cardiovascular diseases, such as atherosclerosis, carotid atheromatous disease, thoracic aortic aneurysms etc., has been discussed by Glenn et al. [22], Bulusu et al. [13,23], and several references cited therein.A common treatment for atherosclerosis is the opening of narrowed arteries resulting from obstructive lesions by angioplasty and stent implantation to restore unrestricted blood flow.The implantation of a realistic "SX-type" stent model as described by Hanus and Zahora [24] led to the observation of transitional secondary structures over those previously reported by Dean [25,26], Lyne [27], Sudo et al. [28], Boiron et al. [29], Timité et al. [30] and Glenn et al. [22].Stent fractures (SFs) may be associated with unanticipated late complications, including clinical in-stent restenosis (ISR), stent thrombosis and aneurysm formation [31][32][33][34][35]. "Type IV" stent fractures have been defined by Jaff et al. [36] as a complete transverse, linear fracture of stent struts along with displacement of the stent fragments.Popma et al. [31] presented excellent angiographic evidence of fracture Types II-IV and noted a high incidence percentage of angiographic restenosis associated with the Type IV fracture (seven out of 18 patients).
Accordingly, the present study is an investigation of stent implantation-related experimental scenarios for the application of entropy-based optimal wavelet scale search and coherent secondary flow characterization, as shown in Figure 1b-d, pertaining to a clean curved artery, a stent-implanted curved artery model and a curved artery with a stent that underwent a "Type IV" fracture.

Objectives
This paper addresses the application of wavelet-decomposition methods, particularly the continuous wavelet transform (CWT) approach in coherent structure identification.
The particle image velocimetry (PIV) technique was used for secondary flows in a model 180-degree curved artery with and without stent implantations to explore the hemodynamics of post facto clinical scenarios.
The wavelet transform is used herein to decompose a two-dimensional vorticity field, ω(x, y), into components at different spatial resolutions, the details of which can be found in Section 3.For the purposes of the analysis wavelet-transformed vorticity, ω = f(x, y) in the sample of such data (X) is considered as a random variable.
Furthermore, the paper also addresses the issue of an optimal wavelet scale search or the optimal scale of the mother wavelet basis, ψ(x, y), that, in a general sense, spans the feature space to adequately represent the signals of interest.Here, the space constitutes vorticity fields generated from PIV experiments on a 180-degree curved artery model and the signals of interest, secondary flow structures (swirling motions).
The flow regimes generally observed in our arterial flow experiments had three complexities associated with them, viz., 1.A 180 • curved tube geometry, guided by arterial networks comprised of varying degrees of curvature, bending, branching and tortuosity; 2. Pulsatility associated with physiological flows in large-scale arteries and; 3. The presence of a stent, giving rise to flow perturbations that interact with the Stokes layer in the flows.
The details of the experiments and the context of the wavelet transform application is presented in Section 2. The wavelet-transformed vorticity signals are then treated with Shannon entropy as a measure to achieve the following two objectives that are described in Sections 3 and 4 (4.1 and 4.2): 1. To implement an entropy-aided search for an optimal wavelet scale; 2. To apply an algorithmic approach in the detection of secondary flow structures in a curved artery model; 3. To create a (non-subjective) threshold-free approach for coherent structure detection in fluid flows.

Description of the Experimental Setup
Pulsatile, physiological hemodynamic flows in a canonical 180-degree curved artery phantom were measured using PIV techniques.The resulting secondary flow morphologies had very interesting spatio-temporal characteristics, such as the breakdown of large-scale secondary flow structures into smaller low vortical strength structures having strain-dominated geometries.
The experimental test section shown in Figure 1 constitutes a rigid, 180 • bent tube and is a canonical representation of arterial curvatures.Two-dimensional flow fields were generated using a two-dimensional particle image velocimetry (PIV) system at several instances of a scaled, carotid artery-based flow rate waveform delivered to the curved artery model.At each instance of time, velocity fields were generated from PIV-based experiments.The vorticity fields were then computed from the two-dimensional velocity field data for further analysis and identification of coherent secondary flow vortical patterns.
The particle image velocimetry setup is shown in Figure 1.These include the inflow conditions, the type of working fluid (a Newtonian blood-analog solution reported by Deutsch et al. [37] and Yousif et al. [38]), bulk flow properties, such as kinematic viscosity, 3.55 ×10 −6 m 2 /s ±2.8% (or 3.55 centistokes) measured at room temperature using a standard Ubbelhode viscometer, the Womersley number (4.2), refractive indices of the acrylic-based test section and the blood-analog fluid and PIV recording parameters and results.
A planar, two-dimensional particle image velocimeter (2C-2D PIV) system was used in the experiment to acquire phase-locked velocity measurements in a 180 • bent rigid pipe model representing an idealized curved artery.The flow rate waveform previously reported by Holdsworth et al. [39] is shown in Figure 2 and has been utilized in the investigation of the three clinically-relevant flow scenarios, as shown in Figure 1.The entire experimental test rig was fabricated using acrylic material (refractive index of ≈1.45) to enable optical access upstream, downstream and within the curved artery test section.The inner diameter of the curved artery model was of 12.7×10 −3 m and was designed to enable optical access at five cross-sectional, planar locations.The laser sheet from the PIV system could effectively illuminate the cross-sectional planes and facilitated the acquisition of secondary flow velocity field data.Acrylic pipes of 12.7 × 10 −3 m inner diameter were attached to both the inlet and outlet of the curved artery model to ensure that fully-developed flow entered the test section.In this study, we are presenting data at the 90 • location in the curved artery model for the physiological inflow conditions shown in Figure 2 for the curved artery, curved artery with stent implantation and curved artery with the Type IV stent fracture (Figure 1b-d).The test section also has provisions for additional measurements at 0-, 45-, 135-and 180-degree cross-sectional planes, as shown in Figure 1.
A viscous, Newtonian, blood-analog fluid was pumped from a reservoir through the acrylic inlet pipe connected to the 180-degree curved test section.In order to minimize the optical distortion of the laser sheet, a working fluid composed of 79% saturated aqueous sodium iodide (921 g per 500 mL of water at 298 K), 20% pure glycerol and 1% water by volume and a refractive index of 1.45 (±3.4%) at 298 K was used.
The blood-analog fluid is circulated in the experimental setup using a programmable pump (Ismatec model BVP-Z) that is controlled by a custom virtual instrument written in LabView through a data acquisition card (National Instruments DAQ Card-6024E).The assumption of treating blood as a Newtonian fluid has been primarily associated with the size of the arterial vessel, as discussed by Zamir [40].When the diameter of the vessel is comparable to the scale of the corpuscular structure of blood and homogenous hematocrit distribution, the Newtonian relations between shear stresses and velocity gradients are untenable in the entire blood flow domain.The vessel diameter in our case (12.7 × 10 −3 m) is much larger than the corpuscular scale, and therefore, the Newtonian assumptions in the core regions of the blood vessels are fairly adequate.However, at the near-wall regions, shear-thinning effects may persist that lead to alterations of the near-wall shear stresses, as shown in our recent publication, van Wyk et al. [41].Since the secondary flow characteristics that we are studying fall fairly well within the core regions of the flow, the Newtonian blood assumptions are valid.
External triggering of the PIV system was provided by the programmable pump instrument control module in order to synchronize measurements at every instance of time in the inflow waveform (Figure 2).The physiological inflow waveform comprises 100 discrete, equispaced samples with a 4-second time period.Up to 200 phase-averaged vorticity fields were generated for each discrete (time) instance on the physiological waveform.Flow field illumination in the 2C-2D PIV system was provided by a dual pulse Nd:YAG laser (wavelength = 532 × 10 −9 m).The recording medium was a CCD camera (LaVision Imager Intense 10 Hz) with a 1376 × 1040 pixel array.Flouro-Max red TM fluorescent polymer microspheres (d p ≈ 7 × 10 −6 m) were used as seeding particles.DaVis 7.2 (LaVision, Goettingen, Germany) software was used for image acquisition and processing.
In order to maintain the similarity between in vivo conditions and laboratory experiments, the inflow waveform was scaled by matching the Reynolds number and the Womersley number (α).The details of waveform-scaling to maintain the in vivo physiological carotid artery Womersley number have been reported previously and are repeated here for clarity [13].The assumed mean diameter of the right and left common carotid arteries (CCAs) is ≈6.3 × 10 −3 m.The physiological inflow waveform from Holdsworth et al. [39] was scaled accounting for the test section tube inner diameter of 12.7 × 10 −3 m, while maintaining the physiological Womersley number (α ≈ 4.2) and maximum Reynolds number (Re max = 1655).It followed that the period of the scaled physiological waveform spanned a time period T = 4 s.The peak flow rate of Holdsworth's physiological waveform (≈23.6 mL/s) and the scaled, new peak flow rate was ≈55 mL/s, as shown in Figure 2 in the scaled physiological inflow waveform.

Stent Geometry and the Idealized Type IV Stent Fracture
In this study, the endovascular stent geometry of the "SX-type" geometry described by Hanus and Zahora [24] was fabricated using a rapid prototyping machine (Objet24 Desktop 3D printer) of 88.9 mm in length and 0.85 mm in strut diameter.This stent geometry (SX-type) can be envisioned as a set of horizontally-and vertically-shifted spirals (springs) knitted from Nitinol wires [13].This embodies a self-expandable stent inserted upstream of the test section (Figure 1c,f).A curved stent with a central angle of 45 • and the same strut diameter (0.85 mm) as the straight stent was also fabricated using the rapid prototyping machine (Figure 1e).The combination of the two stent configurations when implanted in the test section, as shown in Figure 1d, embodied an idealized Type IV fractured stent following the classifications of Jaff et al. [36] and Popma et al. [31].
Our previous publications, Glenn et al. [22], Bulusu and Plesniak [13] and Bulusu et al. [23], can be referenced for further details on our experimental conditions, wherein we reported the PIV post-processing routines, the final window size and overlap, the pulse delays, the number of image pairs acquired toward phase-averaging and statistical convergence, the dynamic spatial range (DSR) and the dynamic velocity range (DVR) that varies with discrete instances of time in the pulsatile inflow waveform.The above information is summarized in Table 1 for succinctness.

Continuous Wavelet Transform Using the Ricker Wavelet
The wavelet transform method utilizes an analyzing function (or mother wavelet) that is localized in space.By convolving a dilated or contracted wavelet with a 2D vorticity field (ω), the wavelet-transformed vorticity (ω) field is generated comprising the coherent structures with a wide range of scales and strengths [12,13].
The application of wavelet decomposition methods to experimental data was demonstrated by Schram and Riethmuller [42], who applied it to digital particle image velocimetry (PIV) data for the vortical structure detection of the leading vortex ring generated in an impulsively starting jet flow.Camussi [3] proposed a method based on the analysis of the local energy content at separated scales, to extract the typical wavenumber associated with structures and, therefore, the typical length-scale of coherent structures.Longo [11] used discrete orthogonal wavelets to study turbulence and intermittency associated with coherent structures in free surface flows (waves) using two-dimensional laser Doppler velocimetry (2D-LDV) measurements.Further, Longo [15] effectively used wavelet analysis to build on the digital PIV measurements of turbulent wave breaking phenomena reported by Melville et al. [43].
A brief description of the implementation of our wavelet-transform approach is presented in this section taking into account the advantages that wavelet decomposition methods offer, such as selectivity in scale and space, retention of the idea of a band-pass filter and an exhaustive vortical pattern search capability.The central idea of our approach is that the PIV-generated measurement ensembles (or images) of vorticity (ω) were treated as the untransformed signals.Wavelet-decomposition methods can be applied to generate transformed signals, i.e., wavelet transformed vorticity fields (ω).
The characterization of complex secondary flow morphologies necessitated the application of a rigorous coherent structure detection method.This led us to develop an advanced wavelet-based coherent structure detection algorithm (PIVlet, v 1.2).
In particular, our implementation involves the use of continuous wavelet transform (CWT), which can be associated with wide, unbounded range of scales and unrestricted characteristics of dilation, translation or rotation [2,3,42,44].It is important to note that Farge et al. [12] were among the first to demonstrate two-dimensional CWT on data from numerical simulations pertaining to near-wall dynamics of three-dimensional turbulent channel flow using a complex-valued function (Morlet wavelet).
We followed Farge et al. [12], who clearly outlined the requirements on the selection and the restrictions placed on the basis function (ψ), summarized in the following points: 1. ψ ∈ L 2 (R n ) and ψ are normalized, usually to preserve the L 2 norm (R = real numbers).2. ψ should be admissible; admissibility implies that the function has zero mean ( ψ should have good localization and smoothness, both in physical and spectral spaces. Our CWT implementation for vortical pattern eductions is intrinsically two-dimensional and employs the two-dimensional Ricker function (Figure 3, Equation ( 1)), which is popularly known as "Mexican hat".The Ricker function was chosen as the wavelet-basis that can be scaled or dilated using a shape parameter ∈ R + , where x ∈ R, y ∈ R. If a wavelet contains enough information about a signal to be represented, the levels of the required resolution reduce in terms of computational complexity.We chose the Ricker wavelet (as represented by Equation ( 1)) for this reason.The two-dimensional Ricker wavelet is considered as smooth in the physical space and resembles the shape of vorticity associated with a vortex [42,44,45].In addition, it bears resemblance to an Oseen vortex, as reported in an earlier study conducted by Schram and Riethmuller [42].Varun et al. [45] pointed out that the Ricker wavelet is a promising kernel that resembles most of the test properties for the Hamel-Oseen solution presented by Schram et al. [46].
The CWT process was initialized with an arbitrary coarse scale, and subsequently, a large number of wavelet scales was generated by changing the scale ( ) to attain finer or coarser scales of the Ricker function (Equation ( 1)).In the current study, a uniform Ricker wavelet was used as a basis function; the family of wavelet functions generated were all isotropic in nature, and the scale factor ( ) dilates the shape of the wavelet initialized by during the CWT process.
An additional parameter ∈ R + was introduced in Equation ( 1) to control the uniformity of the basis function (called the anisotropic factor).In general, the anisotropic factor ( ) initializes the preferential shape of the wavelet before starting the CWT algorithm.The anisotropic factor was set to one in this study to ensure that the 2D Ricker wavelet was isotropic.
CWT of a 2D vorticity field, ω(x, y), can be represented by Equation ( 2), wherein the computation was performed in the Fourier domain using PIVlet 1.2 programmed in MATLAB (MathWork, Natick, MA, USA) for the analysis of PIV data.At each scale, the wavelet captured patterns appropriate for that scale.These scale-specific patterns represent organized coherent structures consisting of widely varying aspect ratios, recognizable smaller-scale structures and sub-structure organizations within large-scale structures.
The application of continuous wavelet transforms offers localization at the cost of computational time that is incurred at wavelet-scale-specific Fourier domain multiplications with the vorticity field.At finer scales, wavelets are sensitive to noise occurring on an otherwise smooth background, leading to a number of redundant coefficients as noticeable distractions from the main results [44].
Therefore, parametric threshold values have to be set such that coherent structures can be detected effectively.The goal of such thresholds is to filter any undesirable structures that form a part of the background noise or contribute to the uncertainty in the flow field.However, the process of setting the threshold is subjective and relies on a priori knowledge of the flow field, which is seldom available.
These challenges led to the question of the optimal scale identification; the question of what wavelet scale when applied to the vorticity field ω(x, y) effectively describes all of the coherent structures of interest in the flow field.We therefore created an autonomous method of optimal scale identification, at which the transformed coefficients best represent the original signal or efficiently and faithfully describe the coherent structures of interest.
The details of the implementation of the CWT algorithm, the idea of Shannon entropy as a measure of the optimum scale of the wavelet basis function and several results of two-dimensional vortical patterns have been reported in Bulusu and Plesniak [13] and Bulusu et al. [23].A brief mathematical discussion is presented in the Section 4, to convey our use of Shannon entropy and its application to data filtering, compression and deconvolution.

Interpretation of Entropy toward Optimal Wavelet Scale Search
Wavelet functions can be dilated into infinitely large or infinitesimally small scales.Therefore, the question of when or how to exit from the process of wavelet transformation naturally arises.Zhuang and Baras [47] suggested that the key to choosing the optimal wavelet basis lies in the appropriate parameterization and an adequate performance measure, in addition to the accurate interpretation of physical phenomena.We incorporated the key ideas of expectation through probability distributions (pdf), as discussed briefly in this section, adopted the notions of Shannon entropy from information theory and blended them into the framework of wavelet-based coherent structure detection.

Threshold-Free Coherent Structure Detection and the Context for the Optimal Wavelet Scale
The context of our interest in the link between wavelet transform and entropy was established by the necessity to resolve swirling (vortical) flow patterns.Accordingly, the objective of this paper is to discuss the application of (Shannon) entropy as a measure in the search of an optimal wavelet scale.Entropy is computed for wavelet scale-specific, wavelet-transformed vorticity and tracked for optimal number of bits, which we hypothesized to sufficiently represent the vortical pattern information in the flow field.
The abstraction of Shannon entropy as a measure of randomness provides a framework to compute the minimum amount of data required to describe any image without loss of information [48,49].Coifman and Wickerhauser [50] derived entropy-based algorithms for selecting the best basis from a library of wavelet packets, guiding our interest in the continuous wavelet transform approach.Coifman and Wickerhauser [50] and their discrete wavelet packet (DWPT) application has been discussed by Wickerhauser [51], Song [52] and implemented by Ruppert-Felsot et al. [53] and Fischer [54].Ruppert-Felsot et al. [53] and Fischer [54] provided a clear description of entropy as a tool to seek a basis in discrete wavelet packets that carries substantial information about the turbulent flow fields that they investigated.
In addition to the aforementioned references, the use of Shannon entropy in this study was drawn from a variety of problems outside of classical information theory, such as image texture characterization and the general problem of data filtering and deconvolution [55][56][57].
Zhuang and Baras [47] have followed Coifman and Wickerhauser [50] and provided a Shannon entropy-based optimal wavelet basis and a gradient-based optimization in discrete wavelet transforms for signal representation.
In the present study, the Shannon entropy of the wavelet-transformed vorticity fields was monitored for each wavelet scale ( ).Zhuang and Baras [47] defined this wavelet-scale-specific entropy as the decomposition entropy.An optimal scale of the basis function is defined as that where coefficients of the transformed signal most efficiently represent the original signal.First, we interpret the use of entropy as a means to capture the amount of randomness or uncertainty in a random variable (X ∈ R), where X can take on any value with its domain, R.
The following brief discussion on the 'expected value' has been incorporated from the book by Cover and Thomas [18] to provide the necessary context in the interpretation of entropy as a measure of uncertainty.The probability of a value of x (trial) is the proportion of the number of times that value is expected over a very large sample (collection of trials, X).Then, p(x) = p(X = x), x ∈ R and p(x) is the probability distribution function (pdf ), such that x∈R p(x) = 1.The statistical description of the random variable data can also be interpreted using the pdf .For example, the notion of the mean of the sample X is described as the expected value: The idea of expected value can be extended to variance and other statistical measures.
In view of the above discussion, we note that Shannon entropy may also be interpreted as a function of the expectation of a discrete random variable I(X), where I(X) = −log 2 p(X).Here, I(X) denotes the information content of X.Hence, the expression for Shannon entropy, H(X), assumes the following form: We set p log 2 p = 0, when p = 0.
For the sake of convenience in the remainder of the paper, logarithms of all quantities are considered to be base two.It should be noted that entropy, H(X) (Equation ( 4)), is a function of the probabilities of the random variable, p(X = x), in a sample (X) and not the random variable itself (x).The entropy associated with each data element of the random variable (x) in the sample (X) is log p(X = x).The expected value of entropy (−E[log p(X)]) of the sample (X) is similar to the expected value (or mean) of the random variable (x) defined in Equation (3), in that we are calculating the expected value (or mean) of the probabilities of the data elements instead of the data itself.Accordingly, the entropy of the sample (X) is the average of the entropies of all of the data elements in the sample, weighted by the probabilities of each corresponding data element.

Algorithm of Optimum Wavelet-Scale Search
Our algorithmic approach relies on the computation of the velocity gradient tensor (matrix in 2D), the generation of vorticity fields (ω(x, y)) and, consequently, wavelet-transformed vorticity fields, ω = f(x, y).We impose the following assumptions: 1. Vorticity data images with greater background noise have greater uncertainty and higher entropy, as they would convey more information than what is necessary to describe vortical patterns.2. Coherent secondary flow structures are associated with clustered pieces of information represented by groups of pixels in vorticity data.3. Coherent structures are associated with certain frequency bands in the energy spectrum of the vorticity fields.4. Background noise is considered to be the flow activity without any measureable or detectable flow characteristics (randomness).
The procedure innately involves the computation of Shannon entropy (Equation ( 4)) from PIV-generated, two-dimensional, wavelet-transformed vorticity fields, ω = f(x, y), where ω is considered to be a random variable in the sample, n.
The focus is then directed to the set of probabilities p i , for each data-entry (or pixel, f i ) in the vorticity field, where Then, the sequence {p n : n = 1, 2, ...} with p n ≥ 0 for all n may be interpreted as a probability distribution function (pdf ) for the sample space n. f i 's correspond to the discrete values of vorticity in this study, and is the normalized square modulus of the i-th element of the signal, with n number of elements.The entropy of this set defined over the pdf is given by Equation (5).
The (decomposition) entropy is computed at each wavelet scale ( ), for the corresponding wavelet-transformed vorticity field, resulting in the following expression: The term ||f i || 2 /||f || 2 is the equivalent measure of probability in the decomposition entropy, as discussed by Zhuang and Baras [47].
Note that f i ∈ R and the chosen wavelet, i.e., the 2D Ricker wavelet, ψ ∈ L 2 (R n ).Recall that the unit of the entropy (H) is expressed in "bits" when the logarithms are of base two [48].Using Shannon entropy as a measure, the number of transformed coefficients of the vorticity field can be monitored for each wavelet scale wherein the number of coefficients in the transformed field is generally a fraction of the original signal (or untransformed vorticity field).The value of e H(f ) is proportional to the number of coefficients necessary to represent the signal for a fixed mean error [47,51,53].The algorithmic representation of the method described in this section is presented here (see Algorithm 1).
The discussion of the results presented in Section 5 is based on the following points: 1.The frequency of the signals under consideration is associated with 2D vorticity fields (ω(x, y)) and represents the spatial frequency characteristics of swirling flow structures.2. Entropy obtains its maximum when the energy of the signal is uniformly distributed in its frequency domain [47].3. The wavelet-transformed vorticity fields (ω(x, y)) are representative of signals with a nonuniform energy distribution in their signal energy spectrum as a whole.The vortical structures in ω(x, y) can be directly associated with certain frequency bands in such an energy spectrum.4. Recall that a lower entropy value, which is attributed to lower signal uncertainty, translates to a higher concentration of the signal energy over the frequency bands associated with vortical structures.Shannon entropy computed at each wavelet scale (i.e., decomposition entropy) takes advantage of the fact of this nonuniform signal energy distribution in ω(x, y) and is representative of signal energy captured at the frequency bands. 5.The process of wavelet transform retains the spatial (2D) frequency bands associated with vortical structures in the energy spectrum of the vorticity fields.
6. Recalling that an optimal wavelet scale was defined as that where coefficients of the transformed signal most efficiently represent the original signal, it can also be attributed to certain frequency bands in which the decomposition entropy is minimized.

Results and Discussion
The approach outlined in Sections 3 and 4 was implemented using a suite of MATLAB programs, collectively called PIVlet (v 1.2).We demonstrate the utility of this approach on a system of numerically-generated Oseen-type vortices for which the analytical solution is known and then proceed toward analyzing experimental data pertaining to a curved artery using three clinical case scenarios shown in Figure 1b-d.

Estimation of Shannon Entropy-Aided Optimal Wavelet Scale on the System of Oseen-Type Vortices
A validation case was generated by constructing a system of elongated Oseen-type vortices of varying sizes.An example of such a system is shown in Figure 4a.These elongated vortices are representative of strain-dominated vorticity fields encountered in curved artery experiments.The canonical Oseen vortex equation is represented by Equation ( 7), where Γ is the circulation, R is the radius of the vortex core and r is the radial coordinate outwardly extending from the core.V θ (r) can also be rewritten in terms of vorticity ( ω), given by Equation (8).At any point with coordinates (i, j), position vector ( r) and velocity, V θ (r) are computed from the vortex cores using Equation (8).
The velocity gradient matrix (D) is then computed for each value of the M × N array in the analytically-generated Oseen-type velocity fields using a first order central difference scheme.The resulting velocity gradient matrices are 2 × 2 matrices and correspond to each entry in the two-dimensional velocity field (Equation ( 9)).

D =
∂Vx ∂x ∂Vx ∂y ∂Vy ∂x ∂Vy ∂y (9) We utilized certain functions in PIVMat 3.01 (Laboratoire FAST., University Paris Sud, University Pierre et Marie Curie, CNRS: France), a MATLAB-based post-processing software, to accurately compute the velocity gradients.The two-dimensional, out-of-plane vorticity (ω z ) was computed from Equation (10) using the appropriate velocity gradients at every data point in the velocity field.
The vorticity "data" in Figure 4a were generated using the aforementioned approach and are also treated as the reference vorticity field.
In Figure 4b,c, we superimposed a Gaussian noise uniformly in two dimensions, such that each figure is representative of 10 dB and 15 dB signal-to-noise ratios (SNR), respectively.Such a noise superimposition generated several smaller scale vortical patterns of randomized strengths and topologies and made coherent structure detection more challenging, i.e., a more stringent test of the proposed technique.We computed the wavelet-transformed vorticity fields (ω) generated by applying continuous wavelet transform using the 2D Ricker wavelet, as described in Section 3, on vorticity fields with 10 dB and 15 dB SNR.We used the Pearson correlation coefficient (R) and peak signal-to-noise ratio (PSNR, dB), in addition to decomposition (Shannon) entropy (H, bits) to facilitate comparisons with the reference vorticity field and the autonomous identification of the optimal wavelet scale ( optimum ).
As a measure of the disparity level between the reference and wavelet-transformed vorticity fields (or images), we computed the Pearson correlation coefficient (R), a statistical measure comparing between two or more image matrices [58,59].Vorticity (ω) and wavelet-transformed vorticity fields (ω) were treated as image matrices of the same order (m × n), and R ω ω was computed using Equation (11).The correlation values are contained in the interval −1 ≤ R ω ω ≤ +1.If two images are identical, the value of R is equal to one; if they are completely uncorrelated or unidentical, R is equal to zero; and if they correlated in a negative sense, i.e., one image is the negative of the other with identical features, R is equal to −1.The Pearson correlation coefficient is computed for the fields represented in Figure 4b,c, for each wavelet scale ( ) with reference to Figure 4a.
The application of wavelet transform leads to scale-wise band-pass filtering.Since the vorticity fields shown in Figure 4b,c) are derived from Figure 4a by superimposing Gaussian noise, the correlation between vorticity fields and wavelet-transformed vorticity fields is an indicator of vorticity data disparity.In the unlikely case of "perfect" filtering of noise, a perfect correlation (R ω ω = 1) between wavelet-transformed vorticity and the reference vorticity data would be the outcome.However, at the optimal wavelet scale, the wavelet-transformed vorticity fields would not only detect the four main Oseen-type vortical structures shown in Figure 4a, but also the smaller scale vortical patterns.The correlation between reference vorticity and the wavelet-transformed vorticity would, therefore, lie in the range, 0 < R ω ω < 1.
The aforementioned optimal wavelet scale was determined by calculating the decomposition (Shannon) entropy (H, bits) of the wavelet-transformed vorticity field at each wavelet scale factor ( ) as described in Sections 4.1 and 4.2.
PSNR (dB) is defined as the ratio between the maximum possible power of a signal and the power of corrupting noise that affects the fidelity of its representation (Equation ( 12)).
In Equation ( 12), Z is the maximum possible pixel value of the image and is related to the bit depth of the encoded image.We converted the vorticity data into equivalent "eight bits per sample" images, and consequently, the value of Z is 255 [55,60].The mean squared error (MSE) as presented in Equation ( 13) represents the cumulative squared error between vorticity (ω) and wavelet-transformed vorticity fields (ω).PSNR represents a measure of the peak error; a low value of MSE indicates a low value error.We computed PSNR at each wavelet scale with the hypothesis that the optimal wavelet scale would correspond to the low values of PSNR.
The aforementioned metrics of analysis, viz.correlation (R ω ω), decomposition (Shannon) entropy (H) and peak signal-to-noise ratio (PSNR), are incorporated in establishing the validity of the optimal wavelet-scale search for vortical pattern detection.Figures 5 and 6 show the results pertaining to Oseen-type vortices with 10 dB SNR and 15 dB SNR Gaussian noise superimposed.
We first discuss Figures 5e and 6e, wherein the decomposition (Shannon) entropy (H) as a function of wavelet scale ( ) is shown.During the initial stages of the wavelet transform, the entropy tends to attain a maximum value, as depicted in both figures.These local maxima were attained due to the effect of wavelet scales that were detecting the energy of Gaussian noise that was uniformly distributed in its corresponding frequency domain.However, as wavelet dilation progresses and higher wavelet scales are realized, the intrinsic nonuniform energy distribution in the wavelet-transformed vorticity fields (ω) is detected.In order to verify this non-uniformity, we computed the decomposition entropy for the Fourier domain ( Ĥ) corresponding to each wavelet scale and plotted them in Figures 5e and 6e.In both of these figures, the optimal wavelet scale, where the decomposition entropy attained the local minimum, was clearly observed.These scales were recorded as optimum = 4.8 for a 10 dB SNR-based vorticity field and optimum = 5.8 for a 15 dB SNR-based vorticity field.These optimal wavelet scales are attributable to certain frequency bands in the corresponding Fourier domain that are representative of the energy associated with vortical structures.
At the optimal wavelet scales ( optimum = 4.8, 5.8), the wavelet-transform vorticity data are shown with sufficient compression to the extent that the four, Oseen-type, large-scale vortical structures were clearly detected.In addition, the ability of wavelet transform to detect surrounding elongated, strain-dominated, smaller scale vortical structures is also observed.Indeed, the smaller scale structures are due to the randomized nature of the Gaussian noise superimposed, which were decomposed to resemble the skewed vortical topologies.These skewed topologies are not physically present within the analytically-generated vortical structures, i.e., they are artifacts.However, the vortical patterns detected in Figures 5b and 6b due to the nature of the noise superimposed are analogous to the turbulence-induced vortical patterns.Furthermore, Figures 5c and 6c show the incoherent structures visualized as (ω − ω) fields that are representative of the background noise that can be considered as flow activity without any measurable or detectable organized flow characteristics.
The output of the algorithmic approach shows that optimal wavelet scales, optimum = 4.8 (Figure 5d, f) and optimum = 5.8 (Figure 6d, f), were achieved and corresponded to correlation values of R ω ω = 0.62, 0.73 and PSNR = 40.09,40.54 dB, respectively.We expected to minimize the PSNR due to the wavelet transform.The results shown in Figures 5f and 6f indicate that in the proximity of the optimal wavelet scales ( optimum = 4.8, 5.8), PSNR values are close to the local minimum.In general, high SNR relates to data with greater signal power or strength.As expected, the low SNR-related vorticity field (10 dB, Figure 4b) has lower correlation (R ω ω = 0.62) with the reference, noise-less vorticity field.In comparison, the other SNR-related vorticity field (15 dB, Figure 4b) has a higher correlation (R ω ω = 0.73).
In summary, the central idea of decomposition entropy ( Ĥ) can now be extended with the algorithmic approach to experimental PIV representing the three clinical scenarios of interest, shown in the schematic drawings, Figure 1b,c,d.

Application of Secondary Flow Structure Detection in the Curved Artery Model
In this section, we present the application of the CWT/entropy approach to PIV data in a 180 • curved artery model.Secondary flow structures form in the curved artery model due to a combination of centrifugal instabilities and adverse pressure gradients.We present the results of secondary flow structure detection at the 90 • location in the model for three clinically-relevant scenarios (see Figure 1b,c,d).
The dilation, translation, breakdown and dissipation of the large-scale secondary flow structures in the curved arteries have the potential to glean clinical insight into cardiovascular hemodynamics, as evidenced by the clinical studies of Sengupta et al. [19], Stonebridge et al. [20] and Stonebridge and Brophy [21].
These secondary flow structures have been referred to deformed Dean-, Lyne-and Wall-type (D-L-W) secondary flow structure morphologies, which have the potential to affect wall shear stresses, which has been known to proliferate cardiovascular diseases, such as atherosclerosis [13,22,23,29].
In Figure 7, the results of the three clinical scenarios are shown, with Figure 7a representing the artery without stent-implantation, Figure 7b depicting the stent-implanted curved artery and Figure 7c the idealized Type IV stent fracture within the curved artery.Since there is no major turbulent source in our flow, a vortex can be intuitively recognized using streamlines.We therefore reported streamlines in our experimental data as a first measure of vortical activity and the presence of large-scale vortices (Figure 7).Four important wavelet scales showing the extrema of the decomposition entropy ( Ĥ) and their corresponding wavelet-transformed vorticity fields (ω) are presented.
In Figure 7a,b,c, we observed trends similar to the system of the analytical Oseen-type vortical system discussed in the earlier Section 5.1.This is evidenced by the existence of local minima in decomposition entropy ( Ĥ) that correspond to the resolution of the large-scale secondary flow structures.
In this paper, we focused on providing an illustrative example of the ability of the threshold-free, autonomous method in detecting vortical patterns at one angular position in the curved artery model (90-degree location) and one instance of time in the physiological inflow waveform (t/T = 0.23) to demonstrate the utility of our algorithm (see Algorithm 1).
Spatio-temporal variations in vortical patterns continuously occur during the pulsatile, physiological waveform.The flow regimes in the 180-degree curved artery model are either laminar or transitional due to upstream stent-induced flow perturbations.The result in Figure 7a,b,c were generated after phase averaging up to 200 realizations to ensure statistical convergence (see Section 2).
The optimal wavelet scale ( optimum) is identified as the wavelet scale ( ) where we observed no change in the morphologies of the large-scale secondary flow structures in the wavelet-transformed vorticity fields (ω) for any further dilation of the wavelet beyond that scale.This phenomenon is observed across all three case scenarios of the curved artery test section.
In the result presented in Figure 7a,b, the persistence of deformed Dean-, Lyne-and Wall-type (D-L-W) vortical structures was clearly observed at the instance of time (t/T = 0.23) during the systolic deceleration.Figure 7a is a baseline case where the stent was not implanted.A triple pair (D-L-W) vortical system is clearly marked in the ω-field corresponding to optimum.Similar observations of the Lyne-type pair were reported by Boiron et al. [29], Timité et al. [30] and Sudo et al. [28] under pulsatile inflow conditions.The third pair of highly-elongated wall-type flow structures ("W") has been reported by Bulusu and Plesniak [13].The D-L-W vortical configuration in the stent-implanted case (Figure 7b) also demonstrated the clear presence of high-strength deformed D-, L-and elongated W-vortices at t/T = 0.23.The W-vortices appeared less elongated and fragmented in comparison to the baseline.Additionally, the movement of D-vortices away from the wall was observed.In the flow scenario under an idealized Type IV stent-fracture condition (Figure 7c), a complete breakdown of D-L-W vortices was observed.These morphologies are completely different from those observed in the baseline (Figure 7a) and stent-implanted case (Figure 7b).Vortical patterns at t/T = 0.23 can be associated with the phenomenon of vortical breakdown, wherein the D-L-W structures exhibit highly dissipative vortical patterns and multi-scale, multi-strength occurrences.

Conclusions
The questions treated by information theory are traditionally related to the areas of data compression and transmission.In this paper, we applied the fundamental concepts of entropy inspired by information theory and explored the uncertainty associated with the signals originating from experimental and synthetic vorticity field data (data source) that were characterized at the recipient end (post-processing).The process of examination by the recipient involved the resolution of data uncertainty and the understanding of the information encoded in the signal.
This paper was based on the hypothesis that coherent secondary flow structures can be encoded over a range of infinitely many wavelet scales.The vorticity data containing interesting morphologies (i.e., swirling flow patterns) were periodically sampled using continuous wavelet transform (CWT).The sampled data were then represented by probability distribution functions, encoded using Shannon entropy.
Vortical patterns in three clinically-relevant flow scenarios under pulsatile, physiological inflow conditions were addressed as multi-scale and multi-strength occurrences that required scale-specific, swirling flow pattern-related characterization.This was made possible using the continuous wavelet transform-based coherent structure detection algorithm, which did not require any user-specified threshold or subjective reasoning.
The two main conclusions of this study are the following: 1.There is an optimal wavelet scale equivalent to the optimal number of bits required to encode coherent structure information in wavelet-transformed vorticity data.This conclusion is supported by an entropy-aided search for an optimal wavelet scale, an algorithmic approach in the detection of secondary flow structures in a curved artery model and led to the creation of an approach for coherent structure detection, which does not require a subjective user-selected threshold.
2. The implantation of prosthetics, such as stents, alters the morphologies of coherent vertical structures, as evidenced by the vortical breakdown phenomenon observed in the idealized Type IV stent fracture case.Such breakdowns have the potential to alter near-wall stress distributions that affect the progression of cardiovascular diseases.The observation was made possible using the Shannon entropy-based wavelet transform method for coherent structure identification presented in this paper.

Figure 1 .
Figure 1.Experimental arrangement and the investigation scenarios (a); particle image velocimetry (PIV) experimental setup with the curved artery test section; (b) clean artery; (c) artery with stent implantation; (d) artery with Type IV stent fracture; (e-f) Schematic drawings of curved and straight stents.

Figure 2 .
Figure 2. Physiological, carotid artery-based inflow waveform reconstructed from the archetypal waveform reported by Holdsworth et al. [39] using 100 discrete values with a 40-ms time interval.

Figure 3 .
Figure 3. Two-dimensional Ricker wavelet at an arbitrary scale.The color gradients are arbitrarily generated by output arguments of the wavelet function (Equation (1)) computed on a uniform X-Y grid.Red and blue contours correspond to regions of high and low output values, which, along with the intermediate output valued-contours, constitute the shape of the wavelet function.

Algorithm 1 :
Optimal wavelet-scale search input : Vorticity field data of dimension m × n with ensemble size N output: Wavelet-transformed vorticity data of dimension m × n with ensemble size N Initialize: A two-dimensional wavelet function of your choice for CWT; Initialize: A predetermined range of scales for wavelet dilation, ∈ R + ; for k ← 1 to N do foreach Vorticity field, N k of the sample N do Construct an m × n matrix -I k ; ; for s ← 1 to do foreach I k of dimension m × n do CWT at wavelet scale, s; Update I k ;

Table 1 .
Experimental and particle image velocimetry (PIV) recording parameters for secondary flows in the 180 • curved tube.