Zimo Wang
email: zimowang@tamu.edu

Faissal Chegdani
email: faissal.chegdani@ensam.eu

Neehar Yalamarti

Behrouz Takabi
email: btakabi@tamu.edu

Bruce Tai
email: btai@tamu.edu

Mohamed El Mansori
email: mohamed.elmansori@ensam.eu

Satish Bukkapatnam
email: satish@tamu.edu

Acoustic Emission Characterization of Natural Fiber Reinforced Plastic Composite Machining Using a Random Forest Machine Learning Model

Keywords: natural fiber reinforced composites, acoustic emission sensor, random forest, smart sensing and in-process monitoring

Natural fiber reinforced plastic (NFRP) composites are eliciting an increased interest across industrial sectors, as they combine a high degree of biodegradability and recyclability with unique structural properties. These materials are machined to create components that meet the dimensional and surface finish tolerance specifications for various industrial applications. The heterogeneous structure of these materials-resulting from different fiber orientations and their complex multiscale structure-introduces a distinct set of material removal mechanisms that inherently vary over time. This structure has an adverse effect on the surface integrity of machined NFRPs. Therefore, a real-time monitoring approach is desirable for timely intervention for quality assurance. Acoustic emission (AE) sensors that capture the elastic waves generated from the plastic deformation and fracture mechanisms have potential to characterize these abrupt variations in the material removal mechanisms. However, the relationship connecting AE waveform patterns with these NFRP material removal mechanisms is not currently understood. This paper reports an experimental investigation into how the time-frequency patterns of AE signals connote the various cutting mechanisms under different cutting speeds and fiber orientations. Extensive orthogonal cutting experiments on unidirectional flax fiber NFRP samples with various fiber orientations were conducted. The experimental setup was instrumented with a multisensor data acquisition system for synchronous collection of AE and vibration signals during NFRP cutting. A random forest machine learning approach was employed to quantitatively relate the AE energy over specific frequency bands to machining conditions and hence the process microdynamics, specifically, the phenomena of fiber fracture and debonding that are peculiar to NFRP machining. Results from this experimental study suggest that the AE energy over these frequency bands can correctly predict the cutting conditions to ∼95% accuracies, as well as the underlying material removal regimes.

Introduction

Natural fiber reinforced plastic (NFRP) composites have recently received considerable interest in the industry [START_REF] Chegdani | Thermal Effects on Tribological Behavior in Machining Natural Fiber Composites[END_REF][START_REF] Chegdani | Multiscale Tribo-Mechanical Analysis of Natural Fiber Composites for Manufacturing Applications[END_REF][START_REF] Holbery | Natural-Fiber-Reinforced Polymer Composites in Automotive Applications[END_REF][START_REF] Al-Oqla | Natural Fiber Reinforced Polymer Composites in Industrial Applications: Feasibility of Date Palm Fibers for Sustainable Automotive Industry[END_REF]. Their high degree of biodegradability and recyclability offers tremendous economic and ecological advantages for sustainable manufacturing. Machining of NFRPs is an essential operation for finishing industrial products from such composite materials. In fact, the viability of fiber-reinforced plastics for various industrial applications is largely determined by the machined surface quality [START_REF] Ramulu | Machining and Surface Integrity of Fibre-Reinforced Plastic Composites[END_REF]. The texture and the finish of the machined surfaces of fiber-reinforced composites have a fundamental bearing on the mechanical properties (strength and hardness) and the geometric tolerance [START_REF] Herrera-Franco | A Study of the Mechanical Properties of Short Natural-Fiber Reinforced Composites[END_REF][START_REF] Ku | A Review on the Tensile Properties of Natural Fiber Reinforced Polymer Composites[END_REF].

While significant investigations of the mechanics of material removal in fiber-reinforced composites have been reported, relatively fewer attempts have been made to characterize the machining of NFRP composites [START_REF] Nassar | Machinability of Natural Fiber Reinforced Composites: A Review[END_REF]. The machining process for fiber-reinforced composites involves complicated material removal mechanisms [START_REF] Chegdani | Mechanics of Material Removal When Cutting Natural Fiber Reinforced Thermoplastic Composites[END_REF] that include a combination of shear deformation of the matrix material, fractures across the natural fibers, fiber pull-out, fiber-matrix detachment through tensile fracture, shearing and sliding along the fiber orientation directions, and delamination of fiber bundles [START_REF] Chegdani | Mechanics of Material Removal When Cutting Natural Fiber Reinforced Thermoplastic Composites[END_REF][START_REF] De Rosa | Acoustic Emission for Monitoring the Mechanical Behaviour of Natural Fibre Composites: A Literature Review[END_REF][START_REF] Herrera-Franco | Mechanical Properties of Continuous Natural Fibre-Reinforced Polymer Composites[END_REF]. Earlier investigations suggest that the cutting mechanisms may highly depend on fiber orientation angle (FOA) and matrix properties [START_REF] Romhány | Tensile Fracture and Failure Behavior of Thermoplastic Starch With Unidirectional and Cross-ply Flax Fiber Reinforcements[END_REF]. Machining of NFRP is further complicated by the multiscale nature of the fibers themselves [START_REF] Chegdani | Multiscale Tribo-Mechanical Analysis of Natural Fiber Composites for Manufacturing Applications[END_REF], as well as the prevalence of significant thermal effects [START_REF] Chegdani | Thermal Effects on Tribological Behavior in Machining Natural Fiber Composites[END_REF].

For example, Fig. 1 shows the machined surface of an NFRP composite with a rough texture resulting from the presence of fiber protrusions spread over the surface even when the machining process was well controlled. These protrusions result whenever natural fibers are cut and/or pulled out across and along their lay direction. As most natural fibers are highly flexible and heterogeneous over multiple geometric scales (e.g., the material composition and mechanical properties vary greatly from macro-to micro-to nano-scales [START_REF] Chegdani | Multiscale Tribo-Mechanical Analysis of Natural Fiber Composites for Manufacturing Applications[END_REF]), the defects on machined surfaces affect the integrity of the finished surface and consequently reduce mechanical properties and geometric and dimensional tolerances, hence impeding the emerging industrial applications of NFRP composites.

Currently, investigations into the machining of NFRP composites rely heavily on offline approaches that employ mechanical property testing and alternative imaging and microscopy instruments to characterize the strength, thermomechanics, surface morphology, and chip formation [START_REF] Babu | Effect of Machining Parameters on Milled Natural Fiber-Reinforced Plastic Composites[END_REF]. While these methods have shed significant light on various modes of material removal and surface modification [START_REF] George | A Review on Interface Modification and Characterization of Natural Fiber Reinforced Plastic Composites[END_REF][START_REF] Valadez-Gonzalez | Effect of Fiber Surface Treatment on the Fiber-Matrix Bond Strength of Natural Fiber Reinforced Composites[END_REF], they cannot fully capture the effects of the enormous statistical heterogeneity and variations in distributions of fibers over multiple scales of the matrix on the material removal mechanisms. In addition, these offline characterization tools are often unwieldy for real-time monitoring of the NFRP machining processes, and a sensorbased in-process monitoring approach is desirable to observe and delineate the diverse material removal mechanisms prevalent during NFRP cutting, to characterize the cutting process in real time, and to provide real-time predictions of the machined surface quality for timely intervention for quality assurance [START_REF] Romhány | Tensile Fracture and Failure Behavior of Thermoplastic Starch With Unidirectional and Cross-ply Flax Fiber Reinforcements[END_REF][START_REF] Teti | Advanced Monitoring of Machining Operations[END_REF][START_REF] Chegdani | Fiber Type Effect on Tribological Behavior When Cutting Natural Fiber Reinforced Plastics[END_REF].

Conventional sensor technologies may have a limited signal-to-noise ratio (SNR) and are often insensitive to the subtle variations in cutting processes at high precision levels, where the cutting dimensions are of the order of a few micrometers and cutting force values are smaller than those in conventional machining processes [START_REF] Lee | Precision Manufacturing Process Monitoring With Acoustic Emission[END_REF]. Acoustic emission (AE), which emanates from the deformation, fracture, and/or friction mechanisms at the various tool, workpiece, and chip interfaces of machining, can be sensitive to monitor and enable us to discern among the various material removal mechanisms in NFRP machining. AE can also provide information related to chip entanglement, crack formation and propagation, and chip breakage in this process [START_REF] Uehara | Identification of Chip Formation Mechanism Through Acoustic Emission Measurements[END_REF][START_REF] Bhuiyan | An Innovative Approach to Monitor the Chip Formation Effect on Tool State Using Acoustic Emission in Turning[END_REF]. As detailed in Section 2, very little, if any, work has been done to study AE monitoring of the NFRP machining process.

Background

Earlier applications of AE for characterizing traditional glass and carbon fiber-reinforced polymer composites employed approaches developed in the nondestructive testing literature for measuring and localizing the damage [START_REF] Margueres | Comparison of Stiffness Measurements and Damage Investigation Techniques for a Fatigued and Post-Impact Fatigued GFRP Composite Obtained by RTM Process[END_REF][START_REF] Davidovitz | Failure Modes and Fracture Mechanisms in Flexure of Kevlar-Epoxy Composites[END_REF][START_REF] Mizutani | Fracture Mechanism Characterization of Cross-ply Carbon-Fiber Composites Using Acoustic Emission Analysis[END_REF]. For machining applications, however, environmental noise/inference needs to be maintained at a level that ensures an adequate signal-to-noise (SNR) ratio. Many of the underlying dynamics of machining processes are complicated, and the signals are often contaminated with a high level of environmental noise. Consequently, how to analyze the nature of the nonstationary AE signals during NFRP machining remains a major challenge. Initial investigations of AE in machining have focused mostly on relating various statistical quantifiers of the signal patterns to predict machining conditions [START_REF] Dornfeld | Acoustic Emission During Orthogonal Metal Cutting[END_REF][START_REF] Kannatey-Asibu | Quantitative Relationships for Acoustic Emission From Orthogonal Metal Cutting[END_REF]. One of the earlier research thrusts was based on applying advanced machine learning approaches to estimate the machining process state using extracted frequency features of AE signals. Kamarthi et al. [26] and Pittner et al. [START_REF] Pittner | Wavelet Networks for Sensor Signal Classification in Flank Wear Assessment[END_REF][START_REF] Pittner | Feature Extraction From Wavelet Coefficients for Pattern Recognition Tasks[END_REF] presented a neural network that extracted features using the wavelet decomposition on a multisensor setup for detecting tool wear. The microdynamics of the AE signal have also been investigated by Bukkapatnam et al. [START_REF] Bukkapatnam | Fractal Estimation of Flank Wear in Turning[END_REF], and a recurrent neural network framework that uses the fractal properties of the attractor of the underlying dynamic system along with other machining parameters as the training inputs was developed to predict flank wear in machining. Furthermore, Rao et al. [START_REF] Rao | Real-Time Identification of Incipient Surface Morphology Variations in Ultraprecision Machining Process[END_REF] applied a recurrent predictor neural network to capture the complex nonlinear and nonstationary evolution of the process underlying the measured vibration and AE signals and detect the incipient surface deterioration in ultraprecision machining processes. Cheng et al. [START_REF] Cheng | Ultra-Precision Machining Process Dynamics and Surface Quality Monitoring[END_REF] applied the Gaussian process with extracted statistical features as well as features based on the nonlinear recurrent analysis for online prediction of the surface morphology. Some statistical analysis approaches [START_REF] Kannatey-Asibu | Quantitative Relationships for Acoustic Emission From Orthogonal Metal Cutting[END_REF][START_REF] Kannatey-Asibu | A Study of Tool Wear Using Statistical Analysis of Metal-Cutting Acoustic Emission[END_REF] have been presented to analyze the AE for monitoring of the cutting tool condition and to relate the AE signal features to the metal cutting process parameters. Lee et al. [START_REF] Lee | Precision Manufacturing Process Monitoring With Acoustic Emission[END_REF] discussed the unique requirements for monitoring precision manufacturing processes and the suitability of AE as a monitoring technique at the precision scale. A statistical modeling approach from our earlier effort [START_REF] Wang | A Dirichlet Process Gaussian State Machine Model for Change Detection in Transient Processes[END_REF] was aimed at handling transient behaviors during ultraprecision machining processes. However, most of these approaches may rely heavily on data-driven modeling and consequently lack explanations for the physical principles connecting the acoustic emission signal to the process microdynamics. Chang and Bukkapatnam [START_REF] Chang | Towards Characterizing the Microdynamics of AE Generation in Machining[END_REF][START_REF] Bukkapatnam | A Statistical Mechanistic Model of Acoustic Emission Generation in Shear Zone of Machining[END_REF] addressed the connection between AE sensor signal waveforms and the microdynamics for machining ductile metals. The model that represents the propagation of AE waveform signals from the plastic deformation in the shear zone suggests that AE energy is highly related to the shear strain and the shear strain rate.

Investigations into the characteristics of AE waveforms from NFRP composite machining, especially those attempting to understand the relationships connecting the measured AE signals and the NRFP machining microdynamics, have received very little attention. This lack of attention is mostly due to the difficulties of characterizing the complicated material removal mechanisms that govern NFRP composite machining. Some research efforts have also investigated the spectral characteristics AE during the machining processes of polymer matrix composite materials [START_REF] Arul | Online Monitoring of Acoustic Emission for Quality Control in Drilling of Polymeric Composites[END_REF][START_REF] Velayudham | Acoustic Emission Based Drill Condition Monitoring During Drilling of Glass/Phenolic Polymeric Composite Using Wavelet Packet Transform[END_REF]. These investigations suggest that AE frequency components are highly related to flank wear and workpiece surface modifications.

The present work investigates the AE characteristics related to the machining processes of NFRP composites using a data-driven analytic approach that not only allows real-time sensing and detection of the variations in the microstructure of NFRPs during machining processes but also provides analytic insights toward explaining the deformation mechanism during machining of such NFRP materials. An experimental study is presented to relate AE signals from NFRP machining to the process conditions and the underlying material removal mechanisms. The experiments were conducted on an orthogonal cutting testbed equipped with a multisensor data acquisition system. To connect the AE signal patterns to the machining mechanisms of NFRP cutting, we employ a random forest machine learning model. This model aims to capture the underlying nonlinear relationship between predictors from AE signals and the estimate of the process condition without sustaining significant overfitting issues. The model was also employed subsequently to enable us to the spectral components of AE signals that are highly sensitive to various material removal mechanisms during the cutting of NFRPs. The remainder of this paper is organized as follows: the experimental setup and the framework of the present approach is presented in Sec. 3; Sec. 4 presents the results of the random forest modeling study and provides insights into the physical origins of the signal components; and Sec. 5 concludes the paper.

3 Experimental and Analysis Approach 3.1 Experimental Setup and Procedure. The experiments were conducted on an orthogonal cutting testbed as shown in the schematic diagram in Fig. 2. The setup consists of two linear sliders and a workpiece holder/clamping vice. The cutting tool setup was attached to one linear actuator (L70, Moog Animatics, Milpitas, CA) through a customized tool holder. A high-torque servomotor drives the linear actuator to ensure a consistent cutting speed. Mainly two types of sensors were mounted on the machine setup: two accelerometers (Kistler Type 8728A500) were mounted on the tool holder and the workpiece holder separately for gathering vibration signals, while an acoustic emission sensor (Physical Acoustics S9225) was mounted on the cutting tool holder for collecting acoustic emission signals during the orthogonal cutting processes. This sensor monitoring system was equipped with a data acquisition system (National Instrument compactDAQ with DAQ module NI-9223) that allows real-time data collection at a 1 MHz sampling rate for each sensor signal. Unidirectional flax fiber-reinforced polypropylene composites (UDF/PP) mounted on the clamp-on vice were used as workpiece samples to investigate the effect of the fiber orientation angle on the cutting process. The UDF/PP samples shown in Fig. 2(c As noted earlier, the purpose of this study is to assess the ability of AE signals to discern various process conditions including changes in the cutting speed and orientation in real time. The ability to discern changes in the process parameters is essential (but not sufficient) for a sensor to be qualified for real-time monitoring applications. This ability is also important in NFRP machining because of the high levels of variation expected in the fiber orientation related to the cutting direction as well as in the speeds during the machining of complex geometries and profiles. Consistent with this notion, the present approach aims to estimate process conditions based on analyzing the spectral characteristics of AE signals gathered from orthogonal cutting experiments using a random forest model, as summarized in Fig. 3. We employed a total of 18 different experimental settings, each formed by a combination of six cutting speeds, v = 2, 4, 6, 8, 10, 12 m/min, and three FOAs, 0 deg, 45 deg, and 90 deg. Experiments were replicated three times at each combination. At the beginning of each experiment, the surface of each NFRPs sample was polished using sand paper with a grit size of ∼15 µm to ensure a consistent initial condition (i.e., smooth surface with minimal deformed protrusions). As the surface was machined, the data acquisition system recorded signals from three channels, viz., two from the vibration sensors and one from an AE sensor, at a sampling rate of 1 MHz. A moving windowed timefrequency analysis was performed to extract frequency features over the specified 125-µs long time windows. To investigate the relationship between the frequency responses of the AE signal and the fundamental cutting mechanisms during the cutting of the NFRPs, we treated the generated time-frequency features as inputs to a random forest model for characterizing the cutting processes with different fiber microstructures. Given the spectral features extracted from the online data during the machining, the trained random forest model was expected to predict the machining conditions, specifically the cutting speed and the FOA in real time. The details of this approach are presented in Secs. 3.2 and 3.3.

Time-Frequency Analysis and Feature Extraction.

Figure 4 shows the representative sensor signals collected synchronously from the AE sensor and the two accelerometers mounted on the tool holder and the workpiece holder, respectively, during one cutting phase. The tool-approaching and the cutting phases are appropriately marked. As the cutting begins, we observe a sudden increase in the amplitude of Accelerometer 2 (attached to the workpiece holder). This increase allows us to isolate the dynamics of the cutting and the noncutting stages, such as tool-approaching and various cutting phases, by distinct changes in the variance of the vibration signal at the onset of cutting.

The AE signal captures the transients, i.e., time-varying frequency patterns associated with the material deformation and/or fractures that occur at irregular intervals during machining [START_REF] Hase | The Relationship Between Acoustic Emission Signals and Cutting Phenomena in Turning Process[END_REF]. Hence, a time-frequency analysis is needed to generate the spectral features of the transient AE signal. Such temporal-spectral features are visually represented using a spectrogram as shown in Fig. 5(a), where the x-axis represents the time index and the y-axis is the frequency range. The energy of various frequency bands is represented using a color map. Such a spectrogram records the frequency components of the signal and their variations over time. The procedures for generating the temporal-spectral features can be summarized as follows:

(1) A sliding window with window width L = 125 is applied to collect a set of AE signals at time index t, i.e., {x t-L+1 , x t-L+2 , …, x t }. (2) Then, the fast Fourier transformation (FFT) is applied to compute the frequency component:

X (t) k = t n=t-L+1 x n e -i2πkn L , k = 0, 1, . . . , L -1 (1)
where {x n } (for

n = t -L + 1, t -L + 2,…, t)
is the time series of the AE signal within sliding window L, and the fast Fourier transformation generates the frequency components

X (t) = [X (t) 1 , X (t) 2 , . . . , X (t) L ] T .
After sequentially generating the frequency component vectors, X (t) 's (for t = L, L + 1, L + 2, …, T), the matrix of the spectrogram can be represented as M TF = [X (L) , X (L+ϑ) , X (L+2ϑ) , . . . , X (T) ], where ϑ denotes the step of the sliding window. Essentially, the spectral energy over the high-frequency components of AE waveforms is related to the underlying microdynamics of the machining processes. The statistical analysis based on a Gini index [START_REF] Menze | A Comparison of Random Forest and Its Gini Importance With Standard Chemometric Methods for the Feature Selection and Classification of Spectral Data[END_REF], as presented in the sequel, suggests that these energy values can serve as important features to discern between various pairs of cutting conditions.

Apart from the spectral energy distribution via the time-frequency analysis, we had considered a few advanced signal analysis approaches, such as wavelet decomposition [START_REF] Bukkapatnam | Local Eigenfunctions Based Suboptimal Wavelet Packet Representation of Contaminated Chaotic Signals[END_REF], Hilbert-Huang transformation (HHT)/empirical mode decomposition (EMD), and intrinsic time-scale decomposition (ITD). However, these methods have specific drawbacks in terms of being able to capture the energies across various bands. For example, as for wavelet decomposition, it has low resolution on high-frequency components, especially those that capture variations in the process microdynamics. A wavelet packet representation could provide fine resolution within specific frequency ranges [START_REF] Bukkapatnam | Analysis of Acoustic Emission Signals in Machining[END_REF]. However, most of the wavelet packet bases (especially, if the orthogonality condition is imposed) have an infinite support at time and/or frequency domain (i.e., every coefficient is estimated by taking all the samples of the AE signals). This complicates the extraction of the energies over various bands. For the HHT/EMD and ITD [START_REF] Iquebal | Change Detection in Complex Dynamical Systems Using Intrinsic Phase and Amplitude Synchronization[END_REF], it is not straightforward to extract frequency components from the intrinsic mode functions (IMFs). These representations introduce severe edge effects when decomposing the AE signals over multiple levels of IMFs, which amounts to losing the data points in the time domain.

Figures 6(a)-6(c) show the spectrograms of the AE and the two vibration signals, respectively, capturing the temporal-spectral features of the tool-approaching and the cutting phases. The spectrograms shown in Fig. 5 were generated by the FFT with a sliding window of 0.125 ms and a 50% time step (i.e., ϑ = 0.5 L), which translates to a temporal resolution of 0.125 ms (x-axis), and we select 200 Hz as the spectral resolution (y-axis). Our choice of the window length and the overlap is based on the following rationale. AE sensor signals are composed of multiple AE pulses released from various events, such as fracture and dislocation avalanche that occur at frequency ranges of over 200 kHz (i.e., these events last just for a fraction of a microsecond up to a few microseconds (see Fig. 10)). We selected a 0.125-ms long time window, which is of sufficient length to achieve the frequency resolution to discern changes in the spectral shape and at the same time, not too long to ignore time variations of AE events (i.e., the transients). Similarly, a 50% was chosen to extract sufficient number of features and at the same time avoiding significant correlation (more specifically, multiple collinearity) among the features to train the random forest model. Taken together, the window length and the overlap allow us to effectively extract the salient time and frequency characteristics of the measured AE signals. The edge effect due to the sliding window may be ignored as the window size (0.125 ms duration) is smaller than durations of the overall recordings (which normally last for seconds) by 4-5 orders of magnitude. The dimensions of temporal-spectral features are illustrated in Note that the multisensor setup allows a partition of the timefrequency features into cutting versus noncutting. Let τ 1 and τ 2 be the start and the end time stamps, respectively, of the segment for the cutting phase (e.g., the start and the end points of the red frame shown in Fig. 6(a)). Here, we denote the cutting phase partition using a series of indicators {δ (t) cut

} = {δ (L) cut , δ (L+ϑ) cut , δ (L+2ϑ) cut , . . . , δ (T) cut }, where δ (t) cut = 1 for t ∈ [τ 1 , τ 2 ]
and δ (t) cut = 0 otherwise. The response vector Y combines the cutting indicator with other control parameters for the experiments, i.e., Y (t) = [δ (t) cut , FOA (t) , v (t) ] T for t = L, L + ϑ, L + 2ϑ . . . , T, where T is the total length of the acquired data, FOA ∈ {0 deg , 45 deg , 90 deg}, and cutting speed v ∈ {2, 4, 6, 8, 10, 12} (units in m/min). Given the time duration T 0 (in seconds) of the recording, here T = T 0 × Fs, where Fs is the sampling rate (Fs = 1 MHz for all experiment recordings). Then, the extracted features, X's, as well as the process responses, Y's, for a single recording Ω can be represented in the following equation:

Ω = [X (t)T , Y (t)T ] T = X (L) X (L+ϑ) X (L+2ϑ) . . . X (T) Y (L) Y (L+2ϑ) Y (L+2ϑ) . . . Y (T) (2)
3.3 Random Forest Classifier for Process Condition Monitoring. The random forest approach [START_REF] Liaw | Classification and Regression by RandomForest[END_REF] is applied to capture empirical relationships between the AE spectral features {X (t) } and {y (t) }, where y (t) is a subset of response {Y (t) } (e.g.,

y (t) = δ (t)
cut for characterizing cutting versus noncutting phases and y (t) = FOA (t) for classifying the fiber orientation). During the training phase, the random forest approach (in Fig. 7) employs samples obtained, with replacement, from the training set {X (t) } and {y (t) } to build multiple tree learning models. The decision rule for each tree learning model was based on the randomly selected features {X (t) } [START_REF] Botcha | Process-Machine Interactions and a Multi-Sensor Fusion Approach to Predict Surface Roughness in Cylindrical Plunge Grinding Process[END_REF][START_REF] Arumugam | Random Forests Are Able to Identify Differences in Clotting Dynamics From Kinetic Models of Thrombin Generation[END_REF]. Unlike conventional regression approaches [START_REF] Afrin | Balanced Random Survival Forests for Extremely Unbalanced, Right Censored Data[END_REF], a random sampling and selection scheme allows a classification model to capture underlying nonlinear relationships without imposing biased structures or overfitting issues. The detailed procedures on training the random forest models are stated as follows [START_REF] Wu | A Comparative Study on Machine Learning Algorithms for Smart Manufacturing: Tool Wear Prediction Using Random Forests[END_REF].

Given a training dataset Ω ′ = {(X (L)T , y (L)T ),(X (L+ϑ)T , y (L+ϑ)T ), . . . , (X (T)T , y (T)T )} (y (t) = FOA (t) for fiber orientation classification and y (t) = δ (t) cut for the cutting versus noncutting phase), the random forest first generates B new training datasets Ω i (for i = 1, …, B), which have the same sample size as the original Ω ′ , with some samples replaced in the dataset (also called bootstrap sampling). Consequently, some observations may be repeated in the new dataset Ω i .

For each bootstrapped sample, the decision tree is constructed by the following procedures: at each node of the tree model, m features are randomly sampled and the best split among these randomly selected predictors is chosen (also called a feature bagging scheme). The splitting criterion at the node is chosen to minimize the sum of the squares, which allows the partition of m predictors into M regions R 1 , R 2 , …, R M . Assume that the splitting variable X j is split at point x j to segment the half plane as R 1 and R 2 , where R 1 ( j, x j ) = {X|X j ≤ x j } and R 2 ( j, x j ) = {X|X j > x j }. The splitting variable X j and the split point x j should be the optimal solution to the following objective [START_REF] Wu | A Comparative Study on Machine Learning Algorithms for Smart Manufacturing: Tool Wear Prediction Using Random Forests[END_REF]:

min j,xj min c1 Xk ∈R1( j,xj) (y k -c 1 ) 2 + min c2 Xk ∈R2( j,xj) (y k -c 2 ) 2 ⎡ ⎣ ⎤ ⎦ (3)
where the constants c 1 and c 2 are the responses in region R 1 and R 2 . Note that the overall expected output/response in each region can be formulated as f (X) = M s=1 c s I(X ∈ R s ). The optimal constant value, c S , is the average of response y (t) in region R s , i.e., c S = average (y (t) |X (t) ). Hence, the constants can be solved by c 1 = average (y (t) |X (t) ∈ R 1 ( j, x j )) and c 2 = average (y (t) |X (t) ∈ R 2 ( j, x j )). The splitting algorithm is stopped once a fully binary tree model of a certain number of levels (in the default setup, 8 levels) is built.

For our case study, the inputs of the feature space are from the extracted spectrogram {X (t) }, and the training outputs {y (t) } are selected from the responses {Y (t) }. For the parameter selection, we chose the number of trees as B = 100 and 25 variables to randomly sample candidates at each decision split for all the cases tested, which results in a minimal out-of-bag error.

The random forest generally performs an implicit feature selection by using a small subset of all variables {X j } for the classification problem. Such subsets are regarded as "strong variables" for superior performance in handling high-dimensional data/features. The associated Gini importance values are the measures that visualize the results from the implicit feature selection of the random forest model and are commonly used as the indicator(s) of the feature relevance to the classification results.

Let p κ = N κ /N be the fraction of the N κ samples from the class FOA = {0 deg , 45 deg , 90 deg} or δ cut = {0, 1} (κ ∈ FOA or δ cut ) of the total ℕ samples at node ξ. The Gini impurity Gini(ξ) is calculated as follows:

Gini(ξ) = 1 - κ p 2 κ (4)
The decrement ΔGini is caused by the splitting and sending of the samples to two subnodes ξ ϖ and ξ ς (with fractions with respect to two subnodes p ϖ = N ϖ /N and p ς = N ς /N). Threshold x j on variable X j is defined as follows:

ΔGini(ξ) = Gini(ξ) -p ϖ Gini(ξ ϖ ) -p ς Gini(ξ ς ) (5)
A search over all variables at node ξ and over all possible segmentation threshold obtains a pair {X j , x j } that allows a maximal ΔGini value. The decrease in Gini impurity is recorded for all node {ξ}'s for all trees B and accumulated for all variable X j 's as:

I Gini (X j ) = B ξ ΔGini {Xj} (ξ, B) (6)
Intuitively, this Gini importance value I Gini indicates how often a particular feature X (t) k was selected for a split and how large its overall discriminative value was for the classification problem [START_REF] Menze | A Comparison of Random Forest and Its Gini Importance With Standard Chemometric Methods for the Feature Selection and Classification of Spectral Data[END_REF]. (t) } and the response {y (t) } = {δ (t) cut }, a random forest model was trained to discern between the cutting and the noncutting phases. Specifically, for both training and testing data sets, 20% of the data were gathered during the cutting operation and the other 80% data during the noncutting operation. The overall data set was split into five subgroups, and each subgroup was tested with the random forest classification. Each data set selected around 1400 samples for training and about 600 samples for testing. We compared the present approach with Fisher's linear discriminant analysis (LDA). At the outset, conventional LDA approaches are effective only when the signal features gathered from every pair of process conditions (e.g., cutting versus noncutting, and FOA) can be separated (in the feature space) using a linear boundary (as in a hyperplane). The LDA models also assume that feature values from different conditions can be grouped into distinct Gaussian clusters with equal variance and that they are well separated from each other. In contrast, a random forest can separate out the signal features from various conditions using complex nonlinear boundaries and does not make any assumptions on the clustering of the feature values. Table 1 summarizes the classification results (in terms of overall accuracy, true positive rate (TPR), and true negative rate (TNR)) from the presented approach versus LDA for the present case study. Here, the TPR measures the proportion of actual samples that are correctly identified as from the cutting condition, and TNR represents the percentage of the samples correctly classified as from the noncutting condition. It is evident that the average accuracy of the random forest model is ∼94%, which is 20% higher than that of LDA. In addition, TPR increases from 0.45 with LDA to 0.84 with the random forest model.

The resulting Gini importance values from the classification problem of the cutting versus noncutting phases are shown in Fig. 8(a). Note that each y-axis label refers to the frequency (kHz) in the temporal-spectral matrix (e.g., the first Gini index relates to the frequency response in 358.4 kHz). A higher Gini importance value of a frequency range indicates a greater contribution of the spectral feature to improve the accuracy. Interestingly, the most significant frequency features determining the accuracy of the random forest based on Gini importance are in the frequency ranges of 350 kHz to 440 kHz. The spectrogram in Fig. 8(b) suggests that more energy is contained in such a frequency range (highlighted in the pink frame) during the cutting compared with the noncutting phase. However, conventional vibration sensors are not suited for monitoring processes at such high frequencies, as they do not have an adequate dynamic response in this range. This high-frequency response highlights the significance of using an AE sensor instead of the more commonly used vibration sensor for fast, high-frequency (>100 kHz) detection of subtle changes during machining at microscales.

Supervisory Monitoring for the Fiber Orientation

Effect. In the next case study, we investigated the effect of fiber orientation on machining microdynamics and the corresponding AE characteristics. Given the extracted time-frequency features {X (t) } and the response {y (t) } = {FOA (t) }, a random forest model was trained to predict the fiber orientation based on the spectral features of AE. Here, spectral values above 100 kHz in every column of the AE spectrogram (Fig. 9) were used as the features because the high-frequency spectral components have the time resolution needed to detect changes in their incipient states. As shown in Fig. 9, the AE time portraits as well as the temporal-spectral features show significant differences between the fiber orientation with 90 deg and the others under the same cutting speed (4 m/ min). In contrast, the signal variance between fiber orientations with 0 and 45 deg are less noticeable in the time domain compared with the difference between their temporal-spectral features plotted in the spectrograms in Fig. 9.

The results of classification accuracy are listed in Table 2 (results are generated based on five computations on the combined data from different cutting speeds). Here, in the confusion matrix, each row Results from previous investigations [START_REF] Herrera-Franco | A Study of the Mechanical Properties of Short Natural-Fiber Reinforced Composites[END_REF][START_REF] Wambua | Natural Fibres: Can They Replace Glass in Fibre Reinforced Plastics?[END_REF][START_REF] Mahdi | A Finite Element Model for the Orthogonal Cutting of Fiber-Reinforced Composite Materials[END_REF][START_REF] Eriksen | Influence From Production Parameters on the Surface Roughness of a Machined Short Fibre Reinforced Thermoplastic[END_REF] into the machinability of fiber-reinforced composites suggest that as the cutting speed increases, the surface integrity drops. This may be due to the viscoelastic property of the composite materials. Also, higher cutting speeds may increase the shear stress on the chip tool interface, which in turn increases the strain/elongation of the fibers nonlinearly. For all these reasons, more fiber failures may occur (fiber breakages, pull-out, and some inefficient shearing causing microfiber failures of the cellulose structure), along with more plastic deformation of the matrix per unit time for higher speed cutting conditions, resulting in more energy released. The evidence can be found in the time portraits of AE signals (see Fig. 10), which were recorded during two experiments with the same FOA but different cutting speeds. Here, the blue lines refer to the AE signal gathered under machining parameters v = 12 m/min and FOA = 90 deg, and the yellow-lined AE signal is from the experiment with v = 2 m/min and FOA = 90 deg. It may be noticed that even under the same fiber orientation angle, a higher cutting speed (12 m/min) may result in a larger AE envelope. That is, when analyzing the cutting mechanisms with respect to different fiber orientations, the cutting speed also needs to be considered.

Clearly, the increase in the AE amplitude due to an increase in the cutting speed may be attributed to the increase in the material shear rates, and hence, sources of AE increase with an increase in the cutting speed.

The performance of the random forest model (in terms of classification accuracy) for identifying the correct FOA using AE signals collected at six different cutting speeds was compared with that of LDA (see Table 3). Three cases of different FOAs (0 deg versus 90 deg, 0 deg versus 45 deg, and 45 deg versus 90 deg) were investigated under each cutting speed. Compared with LDA, the random forest's accuracy for discerning fiber orientations increased by ∼30% (from 65% to 94.8%). As noted earlier, such an improvement is a consequence of the random forest relaxing the strong assumption LDA makes on how features are clustered and distributed. Pertinently, the random forest employs an ensemble of decision trees, such that each tree is trained using a different, This avoids the correlation between selected features. Also, the random forest provides metrics that help with feature selection as well as interpretation of the classification results. Specifically, it generates the quantifiers, such as Gini Index, to statistically indicate the extent to which every feature contributes to the classification accuracy.

v = 2 v = 4 v = 6 v = 8 v = 10 v =
The Gini index values for three cases are listed in Figs. 11(a

Acoustic Emission Signatures and the Fracture Energy

During NFRP Machining. Unlike the brittle failure mechanisms that occur when cutting synthetic fiber (e.g., carbon and glass fiber) reinforced composites with the thermoset polymer matrix, the NFRP sample with the thermoplastic polymer matrix possesses high flexibility and follows ductile deformation behavior when undergoing transverse shear stress during orthogonal cutting [START_REF] Chegdani | Mechanics of Material Removal When Cutting Natural Fiber Reinforced Thermoplastic Composites[END_REF]. For FOA = 0 deg, the cutting process is dominated by the matrix deformation under the shearing zone with some amount of debonding of the fibers along the cutting direction (matrix-fiber shearing) with a few episodes of fiber pull-out and microfiber failures (see Fig. 12). At FOA = 45 deg and 90 deg, however, micro-failures across the fiber axis are prevalent due to the shearing and tensile forces. Also, for the process under FOA = 45 deg, the interfacial shearing along the fiber orientation causes debonding between the fiber and the PP matrix. Moreover, intermittent fracture is one mechanism under FOA = 90 deg that results in fiber tensile failure and/or pull-out toward the cutting direction with a sudden increase in the cutting force.

Figure 13 shows the schematic diagrams of the chip formation of NFRPs under different FOAs. It may be noticed that because of the different orientation of fibers, the microstructure (elementary/ bundle fibers and matrix-fiber bounds) may undergo anisotropic stress/strain, resulting in distinct failure modes of NFRP components during the material removal process. Consequently, the internal failures of the composites may manifest as different modes (e.g., fiber breakage and fiber-matrix debonding). The fracture energy released from such different failure modes may contribute to the generated AE waveform. The following fracture energy analysis of the cutting processes for NFRPs provides a theoretical explanation of the differences in the released fracture energy due to different microstructure variations and cutting speeds under each fiber orientation condition.

As for the orthogonal cutting processes only for thermoplastic composites such as the PP matrix in the tested NFRP samples, the power released during the shear stress can be formulated as follows:

E S = F S • V S b • t u • V = F S b • t u • cos α cos(ϕ -α) (7)
where V S is the shearing speed, V is the cutting speed, b is the chip width, t u is the uncut chip thickness (100 µm), and α is the rake angle. The shear angle ϕ is a function of the rake angle α and the friction angle β, i.e., ϕ = π/4 -1 2 (β -α). Then, the shear force F S around the shearing zone can be formulated as [START_REF] Wambua | Natural Fibres: Can They Replace Glass in Fibre Reinforced Plastics?[END_REF] follows:

F s = A s • τ = t u b sin(ϕ) • τ (8)
Here, τ refers to the shear strength of the matrix. Overall, the energy from the deformation of the thermoplastic composite is

E S = cos α cos(ϕ-α)•sin ϕ • τ.
For different FOAs, the cutting mechanism may vary and in turn the total energy for each FOA may differ. In addition to the energy released from the plastic deformation of the PP matrix due to the combined shear and tensile forces, the following components may contribute to the fracture energy under different FOAs.

For the cutting condition when FOA = 0 deg, the cutting mechanism is related mainly to matrix deforming and very little fibermatrix debonding and/or interfacial shearing of fibers. The matrix shear debonding energy is denoted by W S , with the propagating failure at distance L, given by [START_REF] Kaelble | Theory and Analysis of Fracture Energy in Fiber Reinforced Composites[END_REF]:

W S = 2πr 0 τ 0 L/γ (9)
Here, r 0 is the fiber radius, τ 0 is the maximum boundary shear stress at the fracture, and the stress concentration coefficient γ = (2G/E) 1/2 /r 0 • ln(r 1 /r 0 ), where r `1r 0 is the matrix thickness, G is the matrix shear modulus, and E is the elastic modulus of the fiber.

Based on a square packing geometry assumption, the fracture energy per unit area A of such a condition can be rewritten as [START_REF] Kaelble | Theory and Analysis of Fracture Energy in Fiber Reinforced Composites[END_REF] follows:

W (0) A = W (0) S A = τ 0 E 2G 1 2 L b f 1 (V)f 2 (V) (10)
In contrast, for FOA = 90 deg (shown in Fig. 13(b)), the flax fibers in the workpiece undergo stresses along the fiber-matrix boundaries in the primary shear zone. Hence, the energy released during orthogonal cutting consists of (1) the matrix material deformation in the shear zone, (2) fiber failure (breakage and pull-out from the matrix), (3) interfacial friction between fiber-matrix bonds, (4) chip motion and material entanglement on the rake face, and (5) the interaction between clearance face of the tool and the machined surface. Here we only consider the AE energy released from the sources 1-3. This is because these are considered as primary sources of AE during NFRP machining, and they can be physically related to the underlying microdynamics, especially the material removal mechanisms prevailing under different fiber orientations.

The summation of the friction energy [START_REF] Kaelble | Theory and Analysis of Fracture Energy in Fiber Reinforced Composites[END_REF] between the fiber and the matrix boundary with a constant friction force, τ f = (f/πd)dx, between the fibers and the boundary of the matrix in the effective region/length of l is as follows:

P f = πdτ f l (11)
When a debonded fiber separates off over length l = L f inside the matrix, the friction work W f between the fiber and the matrix can be formulated as follows: The fracture energy for the fiber breakage/ pull-out is represented [START_REF] Kaelble | Theory and Analysis of Fracture Energy in Fiber Reinforced Composites[END_REF] as follows:

W f = - l=0 l=Lf 2πr 0 τ f dl = πr 0 τ f L 2 f ( 12
)
W D = πr 2 o 2E l=LD l=0 σ - 2τ f l r 2 dl (13)
By substituting L = r 0 σ/2τ f , we have

W D = πr 2 o σ 2 L 6E (14)
Considering the cutting speed and the square packing geometry assumption, in addition to the energy released from the matrix deformation, the fracture energy under FOA = 90 deg contains the following: (1) the energy released from fiber tensile failure per unit area A, W D /A = (σ 2 0 /12E) L b f 2 (V) and ( 2) the energy released from friction per unit area A, W f /A = (τ f L f /2r 0 )f 2 (V ). Here, the maximum debonding length L b = r 0 and the volume fraction functions are f

1 (V ) = ln(r 1 /r 0 ) = {ln[(π/V ) 1/2 -1 ]} 1/2 , f 2 (V ) = V, and f 3 (V ) = (V/π) 1/2 [60]. The maximum debonding length is L b = (r 0 /τ f )[(σ b /2) -τ 0 (E/2G) 1/2 f 1 (V )],
and we assume that L f ≈ L b . Note that Table 4 lists the values of all physical properties related to energy released under different types of failure modes. Hence, as the cutting speed v increases and more area A is machined per unit time (A = b • v, where b is the width of the chip and is approximately equivalent to the width of the sample for all experiments). Consequently, more fracture energy is released at higher cutting speeds, resulting in significantly different AE envelopes under different cutting speeds as suggested in Fig. 10.

Relative to connecting the AE responses to various failure modes, earlier efforts [START_REF] Mizutani | Fracture Mechanism Characterization of Cross-ply Carbon-Fiber Composites Using Acoustic Emission Analysis[END_REF][START_REF] Nazmdar Shahri | Damage Evaluation of Composite Materials Using Acoustic Emission Features and Hilbert Transform[END_REF][START_REF] Nimdum | Use of Acoustic Emission to Discriminate Damage Modes in Carbon Fibre Reinforced Epoxy Laminate During Tensile and Buckling Loading[END_REF] 11( f ), the random forest approach is capable of picking up most significant frequency features (within frequency ranges of 260 kHz-280 Hz and 220-240 kHz) to achieve a classification accuracy around 95.3%. Table 5 suggests that these frequency features determining the classification accuracy for FOA = 45 deg versus 90 deg are related to fiber-matrix debonding and fiber pull-out [START_REF] Mizutani | Fracture Mechanism Characterization of Cross-ply Carbon-Fiber Composites Using Acoustic Emission Analysis[END_REF], which are fundamental cutting mechanisms that distinguish the machining processes between FOA = 45 deg and 90 deg [START_REF] Ramulu | Frequency Analysis and Characterization in Orthogonal Cutting of Glass Fiber Reinforced Composites[END_REF]. As for the 45-deg orientation position, the cutting mechanism consists of compression induced interfacial shearing along the fiber orientation direction that causes fiber-matrix debonding; for FOA = 90 deg, the fractures that contribute predominantly to fiber pull-out emerge intermittently (this may be one reason that the AE waveform under 90 deg shows a significant transient pattern compared with the others). Moreover, fiber breakage is prevalent for both FOA = 45 deg and 90 deg cases. This prevalence may be one reason that the AE frequency responses (350 kHz-440 kHz) to the fiber breakage are less significant in discerning the changes between cutting mechanisms under FOA = 45 deg and 90 deg conditions.

Conclusions

This paper presented the framework of an in-process monitoring scheme for machining processes of NFRP materials. By using a random forest model, we showed that AE signals contain relevant information pertaining to the process. The main contributions of this paper are as follows:

• An in-process monitoring scheme with a multiple-sensor setup for gathering vibrations as well as the acoustic emission information during orthogonal cutting of the NFRP has been discussed. Different machining conditions and material removal routes were analyzed using the AE signals, namely, (a) cutting versus no cutting and (b) cutting with different fiber orientations. • A random forest algorithm was applied for supervised monitoring of the processes with accuracy around 95% to distinguish various fiber orientations toward the cutting direction that are related to distinct cutting mechanisms. Results from the random forest approach indicate that the high-frequency information from the captured acoustic emission signals reflects fundamental cutting mechanism changes. • Further analysis shows that the most significant frequency responses determining the accuracy of the random forest model are highly related to the distinct failure modes under various fiber orientation conditions. Hence, the presented approach connects the AE characteristics to the microdynamics of machining NFRPs.

We note that AE features considered in the present work aggregates multiple mechanisms prevalent during the machining of NFRP composites. Process conditions, especially the cutting speed, determine the material deformation rate. Since AE is released from the various material deformation and fracture mechanisms during NFRP machining, the time and frequency characteristics of the AE waveform varies with the cutting speed. The cutting speed also affects the dynamics of fiber deformation and fracture, matrix failure, as well as the actions at the primary and the secondary and tertiary deformation zones. Particularly, this variation in the AE waveform characteristics in the 50-100 kHz band (as noted by the Gini index) can be used to identify the process condition, including the cutting speed and FOA. The AE waveform, as gathered, is a composition of elementary AE pulses released from various different sources. Toward using the results of these predictions for process control, it would be desirable to isolate contributions of individual source in the measured AE signal. Also, for future work, other factors, such as the effect of fiber (bundle) distributions and chip formations, need to be considered for further AE analysis. Further advances in the analytical approach is needed for identifying how various transient components from specific sources are composed in the measured AE signals to diagnose incipient anomalies during the machining of NFRPs. We also note that the present approach can be applied for predicting process conditions during machining of different types of composite materials. Further experimental studies within machining of alternative composites materials (e.g., synthetic fiber-reinforced plastics or biocomposites) are needed to identify the chief material removal and surface modification mechanism, characterize the resulting AE waveforms, as well as to develop practical methods to monitor in machining process performance and surface quality in real time.

Fig. 1 (

 1 Fig. 1 (a) Surface obtained from polishing a flax fiber-reinforced plastic material showing a clear contour of flax fibers and basin PP matrix (with a magnified SEM image for showing the elementary fiber and the fiber bundle) and (b) surface finished using orthogonal cutting processes in which the torn/broken-off fibers are randomly oriented on the surface finish (with the SEM image showing detailed torn fibers on the finished surface)

  ) are manufactured by the material supplier Composites Evolution-UK. Each UDF/PP workpiece is a rectangular slab with dimensions of 20 × 15 × 4 mm 3 . The fiber volume fraction is 40%, and the fiber unidirectionality is maintained by polyester weft fiber with low-volume fraction (around 5%). The orthogonal cutting experiments were conducted with the depth of cut around 100 µm in different cutting conditions with varying cutting speeds (v = 2, 4, 6, 8, 10, and 12 m/min) and different FOAs (FOA = 0 deg, 45 deg, and 90 deg) relative to the cutting direction (shown in Fig. 2(d)). Polycrystalline diamond (PCD) cutting tool inserts (Sandvik Coromant-Model TCMW16T304FLP-CD10) with a tungsten carbide substrate (as shown in Fig. 2(c)) were employed for the experiments.

Fig. 2

 2 Fig. 2 Experimental setup for the orthogonal cutting process: (a) schematic diagram showing the mounting locations of AE and vibration sensors, (b) the orthogonal cutting machine setup with the installed high-speed camera and illumination system, (c) a PCD inserted cutting tool is held against the unidirectional flax fiber/PP matrix sample before the cutting begins, and (d) a schematic diagram showing the disposition of the cutting direction relative to the fiber orientation angle

Fig. 4

 4 Fig. 4 Time portraits showing the synchronously gathered AE signal along with vibration signals from the tool holder (Accelerometer 1) and the workpiece holder (Accelerometer 2) during an orthogonal cutting process experiment (v = 8 m/min, FOA as 45 deg)

Fig. 5 .

 5 The energy of each element in the spectrogram matrix is represented in a color map, ranging from -70 dB (blue colored) to 0 dB (red colored). The sensitivity of the AE sensor can be validated by the high-frequency response of the AE signal (Fig. 6(a)) compared with the vibration signals shown in Figs. 6(b) and 6(c). In fact, only the spectrogram of the AE signals showed remarkable/noticeable differences in the high-frequency ranges (>100 kHz).

Fig. 5 A

 5 Fig. 5 A schematic diagram showing the temporal-spectral feature extraction with time stamp resolution of 0.125 ms, a 0.6725 ms overlap, and the frequency resolution as 200 Hz. The process responsesY (t) = [δ (t)cut , FOA (t) , v (t) ] T , where FOA = 0 deg and v = 4 m/min. Note that the matrix M TF records the magnitudes of the Fourier components (as Volts) and their logarithmic amplitudes (in decibels) are depicted in the spectrogram to visually signify the frequency components' variations.

Fig. 6

 6 Fig. 6 The representative time-frequency domain spectrograms (in dB) gathered during an orthogonal cutting experiment (v = 8 m/min, FOA as 45 deg): (a) the AE sensor, (b) the Accelerometer-1, and (c) the Accelerometer-2. Here, the rectangular (red) frame at the right side of every plot indicates the temporal-spectral information during the cutting phase.

Fig. 7

 7 Fig.7The schematic diagram of the random forest machine learning approach

Fig. 8

 8 Fig. 8 Results for classifying cutting versus noncutting: (a) the Gini importance plot with the significant features labeled (y-axis) by the frequency band (kHz) (e.g., the highest Gini value corresponds to the frequency response at 358.4 kHz) and (b) a spectrogram representation of the noncutting phase versus the cutting phase (white line segmented) with corresponding significant frequency bands highlighted in the pink frame

Fig. 9

 9 Fig. 9 The combined AE signals recorded during orthogonal cutting with different FOAs at the same cutting speed (4 m/min) and three spectrograms, one for each FOA

Fig. 10

 10 Fig. 10 Comparison of the waveforms of the AE signals between two cutting speeds: 2 m/min (yellow) versus 12 m/min (blue) for FOA = 90 deg

  )-11(c); Figs. 11(d)-11( f ) depict comparisons of the spectrograms corresponding to these three cases (0 deg versus 45 deg, 0 deg versus 90 deg, and 45 deg versus 90 deg), where the vertical red line in each figure segments two orientation conditions. Note that here the pink blocks highlight the most significant frequency features (with highest Gini index values) suggested by the results from the random forest classification. Early work [51,52] suggests that AE signals with a frequency response above 350 kHz during failure tests of fiber-reinforced composites are more related to fiber breakage. Our case clearly shows that the most significant spectral features that distinguish 0 deg from 45 deg and 90 deg are above 350 kHz. As in the power spectrums of both Figs. 11(d ) and 11(e), more energy in high frequencies (above 350 kHz) is observed for FOAs = 45 deg and 90 deg compared to FOA = 0 deg. When we compare the cutting processes for 45 deg versus 90 deg, the differences in the frequency components are less discernible, but the random forest model is able to select significant frequency features (frequency bands 260 kHz-280 Hz and 220-240 kHz, as highlighted in Fig. 11( f )) to recognize different fiber orientations (45 deg versus 90 deg) for cutting NFRPs. Such different features in certain high-frequency bands suggest that different responses of the AE sensor may be related to different types of failure modes/cutting mechanisms.

Fig. 11

 11 Fig. 11 GINI indices of the significant frequency components for three FOA classifiers: (a) 0 deg versus 45 deg, (b) 0 deg versus 90 deg, and (c) 45 deg versus 90 deg. The spectrograms of the juxtaposed AE signals show the differences in the frequency contents among these three cases: (d) 0 deg versus 45 deg, (e) 0 deg versus 90 deg, and (f ) 45 deg versus 90 deg. The most significant frequency features are highlighted by pink frames.

Fig. 12 A

 12 Fig.12A schematic diagram showing the cutting mechanisms under different fiber orientation angles[START_REF] Ramulu | Frequency Analysis and Characterization in Orthogonal Cutting of Glass Fiber Reinforced Composites[END_REF]

  Random Forest Model to Detect the Start and End of the Cutting Process Using Acoustic Emission Signals. Given the extracted AE time-frequency features {X

	4 Acoustic Emission-Based Monitoring of the NFRP
	Machining Process
	4.1 A

Table 1

 1 Cross-validation comparisons of random forest versus LDA based for the study case of cutting versus noncutting

				True positive	True negative
			Accuracy	rate	rate
	Random	Mean	0.9427	0.8390	0.9819
	forest	Std	0.0124	0.0455	0.0079
	Fisher's LDA Mean	0.7674	0.4443	0.8189
		Std	0.0264	0.0566	0.0340

Table 2

 2 Cross-validation results for various fiber orientations using the random forest

				Predicted
	Confusion matrix		0 deg	45 deg	90 deg
	Actual	0 deg	0.8531	0.1259	0.0210
	45 deg	0.1354	0.7446	0.1200
	90 deg	0.2047	0.1913	0.6040

Table 3

 3

	Accuracy of FOA classification for various fiber
	orientations under each cutting speed (m/min) (in comparisons
	with LDA)

Table 5 .

 5 included extensive mechanical tests on the composite materials and applied frequency analysis to the resulting acoustic emission signals for different types of failures. A summary of the frequency responses of AE signals to different failure mechanisms of fiber-reinforced composites is listed in Matching significant frequency features (results from random forest classifications with the highest Gini coefficients) with the AE characteristics (frequency responses) for different failures modes allows us to identify different failure modes related to the prevalent cutting mechanisms under distinct fiber orientations. As shown in Figs. 11(a), 11(b), and 11(d), the AE frequency responses ranging from 350 kHz to 440 kHz distinguish FOA = 0 deg from other angles. That is, compared with the fiber debonding (which has an AE frequency response in the range from 240 kHz to 310 kHz), which is prevalent under orthogonal cutting with 0 deg FOA, the frequency responses related to fiber breakage (350 kHz to 440 kHz) are more significant for both FOA = 45 deg and 90 deg. When we compare the results of FOA = 45 deg versus 90 deg, even though differences in frequency components are less noticeable as shown in Fig.

Table 4

 4 Physical properties of the NFRPs materials ≈ 0.7676, f 2 (V ) = 0.4, and f 3 (V ) = 0.3568 Table5Dominant range of frequency response (kHz) to fracture energy[START_REF] Nimdum | Use of Acoustic Emission to Discriminate Damage Modes in Carbon Fibre Reinforced Epoxy Laminate During Tensile and Buckling Loading[END_REF][START_REF] De Groot | Real-Time Frequency Determination of Acoustic Emission for Different Fracture Mechanisms in Carbon/Epoxy Composites[END_REF][START_REF] Park | Interfacial Evaluation and Durability of Modified Jute Fibers/ Polypropylene (PP) Composites Using Micromechanical Test and Acoustic Emission[END_REF]

	List of physical properties
	Fiber radius r 0	∼10 µm
	Fiber tensile strength σ 0	∼500 MPa [54,55]
	PP matrix tensile strength	∼29.5 MPa [55,56]
	Maximal shear stress debonding	∼2 MPa [57]
	Fiber Young modulus E	∼50 GPa [55]
	PP matrix Young modulus E	∼1 GPa
	Constant frictional stress τ f	7-40 MPa
	Possible range of shearing stress	30 MPa [58]
	Fiber volume fracture V	0.4
	Fiber packing	Square (4 nearest neighborhoods)
	Volume fraction functions	f 1 (V)

Acknowledgment

The authors thank the National Science Foundation (CMMI #1432914) for supporting the research.


Clement Papadacci

Victor Finel

Jean Provost

Olivier Villemain

Patrick Bruneval

Jean- Luc Gennisson

Mickael Tanter

Mathias Fink

Mathieu Pernot
email: mathieu.pernot@inserm.fr

Imaging the dynamics of cardiac fiber orientation in vivo using 3D Ultrasound Backscatter Tensor Imaging

Keywords: ultrasound, ultrafast imaging, myocardium, echocardiography, anisotropy

The assessment of myocardial fiber disarray is of major interest for the study of the progression of myocardial disease. However, time-resolved imaging of the myocardial structure remains unavailable in clinical practice. In this study, we introduce 3D Backscatter Tensor Imaging (3D-BTI), an entirely novel ultrasound-based imaging technique that can map the myocardial fibers orientation and its dynamics with a temporal resolution of 10ms during a single cardiac cycle, non-invasively and in vivo in entire volumes. 3D-BTI is based on ultrafast volumetric ultrasound acquisitions, which are used to quantify the spatial coherence of backscattered echoes at each point of the volume. The capability of 3D-BTI to map the fibers orientation was evaluated in vitro in 5 myocardial samples. The helicoidal transmural variation of fiber angles was in good agreement with the one obtained by histological analysis. 3D-BTI was then performed to map the fiber orientation dynamics in vivo in the beating heart of an open-chest sheep at a volume rate of 90 volumes/s. Finally, the clinical feasibility of 3D-BTI was shown on a healthy volunteer. These initial results indicate that 3D-BTI could become a fully noninvasive technique to assess myocardial disarray at the bedside of patients.

Introduction

The myocardial fibers architecture plays a major role in the cardiac function. Fibers orientations vary continuously and smoothly through myocardial walls [START_REF] Streeter | An engineering analysis of myocardial fiber orientation in pig's left ventricle in systole[END_REF] and this complex organization is closely linked to the mechanical and electrical myocardial function [START_REF] Arts | Relating myocardial laminar architecture to shear strain and muscle fiber orientation[END_REF][START_REF] Costa | Laminar fiber architecture and three-dimensional systolic mechanics in canine ventricular myocardium[END_REF][START_REF] Waldman | Relation between transmural deformation and local myofiber direction in canine left ventricle[END_REF] . For instance, the helicoidal fiber distribution in the left ventricular walls contributes to the torsional motion of the heart during the ejection phase, which is believed to improve the heart's pumping efficiency [START_REF] Taber | Mechanics of ventricular torsion[END_REF][START_REF] Hansen | Effect of volume loading, pressure loading, and inotropic stimulation on left ventricular torsion in humans[END_REF] . The fibers of the heart are also related to its electromechanical properties [START_REF] Hooks | Laminar arrangement of ventricular myocytes influences electrical behavior of the heart[END_REF][START_REF] Kadish | Interaction of fiber orientation and direction of impulse propagation with anatomic barriers in anisotropic canine myocardium[END_REF][START_REF] Roberts | Influence of cardiac fiber orientation on wavefront voltage, conduction velocity, and tissue resistivity in the dog[END_REF][START_REF] Taccardi | Effect of myocardial fiber direction on epicardial potentials[END_REF] ; indeed, the electrical activation propagates preferentially in the direction of myocardial fibers and results in their synchronous contraction. Myocardial fiber disarray is thought to appear in the early stage of many pathologies such as in cardiomyopathies or in fibrosis [START_REF] Tseng | Imaging myocardial fiber disarray and intramural strain hypokinesis in hypertrophic cardiomyopathy with MRI[END_REF] . The mapping of the fiber architecture, on the one hand, could increase our knowledge of the cardiac function [START_REF] Savadjiev | Heart wall myofibers are arranged in minimal surfaces to optimize organ function[END_REF] , and on the other hand, could potentially enable early diagnosis of cardiomyopathies.

Yet, no method to map the fiber architecture is currently used on a regular basis for clinical purposes.

Magnetic Resonance Diffusion Tensor Imaging (MR-DTI) [START_REF] Hsu | Magnetic resonance myocardial fiber-orientation mapping with direct histological correlation[END_REF][START_REF] Scollan | Histological validation of myocardial microstructure obtained from diffusion tensor magnetic resonance imaging[END_REF] is widely used for mapping the structural connectivity of the human brain. However, its application to the heart, while feasible in animals [START_REF] Reese | Imaging myocardial fiber architecture in vivo with magnetic resonance[END_REF] and humans [START_REF] Wu | Diffusion Tensor Magnetic Resonance Imaging Mapping the Fiber Architecture Remodeling in Human Myocardium After Infarction: Correlation With Viability and Wall Motion[END_REF] , remains challenging in vivo due to limited frame rates (and hence long acquisition times) and limited robustness of the MR-DTI signals to tissue motion [START_REF] Tseng | Imaging myocardial fiber disarray and intramural strain hypokinesis in hypertrophic cardiomyopathy with MRI[END_REF][START_REF] Reese | Imaging myocardial fiber architecture in vivo with magnetic resonance[END_REF][START_REF] Wu | Diffusion Tensor Magnetic Resonance Imaging Mapping the Fiber Architecture Remodeling in Human Myocardium After Infarction: Correlation With Viability and Wall Motion[END_REF] .

Optical methods such as optical coherence tomography [START_REF] Fleming | Quantification of cardiac fiber orientation using optical coherence tomography[END_REF] and two-photon microtomy [START_REF] Huang | Three-dimensional cardiac architecture determined by two-photon microtomy[END_REF] can also map the fiber directions at the microscopic level but remain limited to superficial ex vivo tissue and to small regions of interest. In ultrasound imaging, several methods have been proposed to quantify the anisotropy of various physical parameters linked to the fiber orientation, including the ultrasonic attenuation [START_REF] Mottley | Anisotropy of the ultrasonic attenuation in soft tissues: measurements in vitro[END_REF][START_REF] Baldwin | Measurements of the anisotropy of ultrasonic attenuation in freshly excised myocardium[END_REF] , the integrated backscattered intensity [START_REF] Mottley | Anisotropy of the ultrasonic attenuation in soft tissues: measurements in vitro[END_REF][START_REF] Baldwin | Measurements of the anisotropy of ultrasonic attenuation in freshly excised myocardium[END_REF][START_REF] Wickline | Three-dimensional characterization of human ventricular myofiber architecture by ultrasonic backscatter[END_REF][START_REF] Madaras | Anisotropy of the ultrasonic backscatter of myocardial tissue: II. Measurements in vivo[END_REF] and the myocardial stiffness using Elastic Tensor Imaging (ETI) [START_REF] Lee | Mapping myocardial fiber orientation using echocardiography-based shear wave imaging[END_REF] . However, none of these methods has been yet implemented for the time-resolved, 3D mapping of the myocardial fibers orientation.

In this paper, an entirely novel ultrasound-based imaging technique called Backscatter Tensor

Imaging (3D-BTI) is described and evaluated in vivo with respect to its capability of mapping the myocardial fibers orientation and its time-resolved dynamics during an entire cardiac cycle. 3D-BTI is based on ultrafast volumetric ultrasound acquisitions, which are used to quantify the spatial coherence of backscattered echoes at each point of the imaged volume at a rate of 90 volumes/s. Spatial coherence is a physical property of backscattered waves that depends on the distribution of scatterers at the subwavelength scale. Ultrasonic spatial coherence has thus the capability to provide statistical information on the tissue microstructure at a scale far below the resolution of ultrasound images. In 3D-BTI, the anisotropy of spatial coherence is analyzed on a 2D matrix probe in order to derive its principal directions and derive the fibers orientation.

The link between the anisotropy of spatial coherence and the fibers orientation was previously demonstrated in composite materials and in fibrous soft tissues such as the myocardium and the skeletal muscle [START_REF] Derode | Spatial coherence of ultrasonic speckle in composites[END_REF][START_REF] Papadacci | Ultrasound backscatter tensor imaging (BTI): analysis of the spatial coherence of ultrasonic speckle in anisotropic soft tissues[END_REF] .

Another key aspect of 3D-BTI is the use of volumetric ultrafast plane wave imaging. We recently developed a 3D fully programmable ultrasound system for volumetric ultrafast plane wave imaging achieving volume rates of up to 5000 volumes/s [START_REF] Provost | 3D ultrafast ultrasound imaging in vivo[END_REF] . By coherent compounding of the backscattered echoes associated with tilted plane wave emissions, the spatial coherence and thus the fibers orientation can be obtained simultaneously in each voxel of the volume at high volume rate. 3D-BTI is hereby found to be capable of mapping the transmural fibers orientation in a beating heart. For the purpose of this validation study, 3D-BTI was compared to histological analysis on ex vivo explanted myocardial tissues. An excellent agreement on the transmural fiber orientation was found between 3D-BTI and histology. The in vivo feasibility of 3D-BTI on a beating heart was then demonstrated in an open chest ovine model. Finally, transthoracic imaging was performed on a healthy volunteer to show the clinical feasibility of 3D-BTI. To our knowledge, it is the first time that the dynamics of the myocardial fibers are observed noninvasively over an entire cardiac cycle of a single heartbeat in humans.

Results

The general principle of 3D-BTI is illustrated in figure 1. First, tilted plane waves (figure 1 A) are emitted from a 2D matrix array probe connected to a customized, programmable, ultrasound system [START_REF] Provost | 3D ultrafast ultrasound imaging in vivo[END_REF] . For each emitted plane wave, the backscattered echoes received by each element of the matrix probe are recorded and further processed using 3D plane-wave coherent compounding to synthetically generate a voxel-specific focal region (figure 1 B). The spatial coherence associated to each focal region is then computed at each point of the 3D volume (figure 1 C) and used to determine the fiber orientation by applying an elliptic fit (figure 1 D).

Finally, a vector representation is used for the visualization of fibers in 3D space (i.e. figure 1 E).

Ex vivo experiments

To establish the accuracy and precision of 3D-BTI, myocardial samples explanted from the anterior wall of the left ventricle of a porcine heart were imaged using 3D-BTI and further analyzed by histology. figure 2.A shows the 3D grayscale ultrasound image of one myocardial sample and three examples of spatial coherence functions obtained at different depths (figure 2 B, C, D). At each location, the spatial coherence presents a strong anisotropy with a principal direction that varies with depth. The vector representation (figure 2.F) of the entire volume shows that fiber orientation varies across the wall. Fiber orientation was found to vary gradually through the ventricular wall, which is consistent with the literature 1 . The transmural orientation was averaged over the volume for the 5 myocardial samples and was found to vary continuously through the wall with an average difference of 98.6°± 8.9° between endocardium and epicardium. These results were compared against the fiber orientation obtained from histological analysis.

In vivo open chest experiments

The in vivo feasibility of 3D-BTI was demonstrated in the beating heart of an open chest sheep.

3D-BTI acquisitions were performed during a complete single cardiac cycle with a volume rate of 90 volumes/s. Figure 4 shows examples of fiber reconstruction at specific times of the cardiac cycle (A Late Diastole, B Early Systole, C Late Systole, D Early Diastole). It demonstrates the feasibility of following the fiber orientation during an entire cardiac cycle at high frame rate (the complete cine-loop is shown in supplemental movie 1). The transmural fiber distribution at one location as a function of time is shown in figure 5. Fiber angles were found to vary as a function of depth with an absolute difference of 96° from the epicardium to the endocardium.

In vivo transthoracic BTI of the human heart

Finally, the feasibility of 3D-BTI on the human heart was shown on a healthy volunteer in a realistic clinical transthoracic imaging setting. Transthoracic conventional gray-scale imaging was used to position the probe on the parasternal view. 3D-BTI was then performed with a total imaging depth of 60 mm. Figure 6 shows the helicoidal fiber distribution imaged in the left ventricle anterior wall at mid-level both in systole and in diastole. Fiber angles were found to vary transmurally with an absolute difference of 104±9° from the epicardium to the endocardium.

figure 6 Transthoracic imaging of myocardial fiber orientation in the human heart.

The fiber orientations of the antero-septal LV wall is shown in end-diastole and end-systole.

Discussion

In this study, we have introduced 3D-BTI, a novel ultrasound-based technique for the mapping of the fibers orientation in myocardial tissues. The objectives were to determine the accuracy and precision of 3D-BTI in mapping the fibers structure and to demonstrate its feasibility on beating hearts. The accuracy and precision of 3D-BTI were validated in 5 ex vivo porcine left ventricular myocardial samples. The fibers angles distribution of each sample was obtained within volumes of 3 cm 3 and the transmural angle variations of each sample were compared to histology. An excellent agreement was found against histology, therefore demonstrating the accuracy and precision of 3D-BTI.

The in vivo feasibility of 3D-BTI was then demonstrated in a beating sheep heart. The fibers distribution of the left ventricular free wall was successfully imaged during one entire cardiac cycle at a volume rate of 90 volumes/s. To our knowledge, this is the first time that the dynamics of the myocardial fibers were mapped during an entire cardiac cycle at high frame rate. The transmural fiber variation was found to be 96° in diastole and did not change significantly during the systolic phase despite a 56% increase in wall thickness. This small variation over the cardiac cycle confirms the results of several studies that have investigated the fibers orientation in diastole and in systole using histology [START_REF] Streeter | Fiber orientation in the canine left ventricle during diastole and systole[END_REF] . Finally, the feasibility of BTI in a realistic clinical transthoracic imaging setup was shown on a healthy volunteer and allowed for the mapping of the myocardial fibers orientation of the human heart both in systole and in diastole.

3D-BTI could not only bring new insights in the knowledge of the cardiac mechanics but it could also become a major tool to non-invasively detect fiber disorders linked to cardiac diseases in early stages. Myocardial fiber disarray has been shown in post-infarct remodeling of the left ventricle as well as in hypertrophic cardiomyopathy. 3D-BTI may be used to quantify disarray in the fiber angles but also to map the fractional anisotropy in order to quantify collagen infiltration and fibrosis content. Fractional anisotropy provides a measure of the level of anisotropy and can be obtained using the coherence values along and across the fibers.

3D-BTI presents a number of advantages as it can be performed non-invasively at high frame rate, and in vivo in beating hearts. In contrast to shear wave imaging based techniques, 3D-BTI is based on low energy emissions and thus can be performed continuously to provide the complete description of the fiber dynamics at high frame rate. Moreover, while the present study focused on the use of plane waves for the generation of synthetic foci, diverging waves [START_REF] Papadacci | High-contrast ultrafast imaging of the heart[END_REF][START_REF] Provost | A clinical feasibility study of atrial and ventricular electromechanical wave imaging. Heart Rhythm Off[END_REF] , which can achieve even larger fields of view and image the entire heart in 3D [START_REF] Provost | 3D ultrafast ultrasound imaging in vivo[END_REF] , could also be used to improve 3D-BTI. The use of diverging waves would however result in a non-uniform spatial resolution so that this limitation needs to be evaluated in more details.

3D-BTI is robust to motion as a result of the use of ultrafast volumetric acquisitions, which enabled the mapping of the coherence functions in entire 3D volumes using a limited number of transmits. For example, while 49 plane wave emissions were sufficient to perform 3D-BTI in each voxel, performing the equivalent acquisition using standard focusing would have required 64*64*500 emissions, which would corresponds to a few minutes per volume instead of a few milliseconds in the case of plane wave imaging.

One of the current limitations of 3D-BTI is its incapability to map the z component of the fiber direction. Despite this limitation, mapping of fiber directions remains possible in anisotropic soft tissues in which the fibers are parallel to the 2D array plane, which is the case for the skeletal muscle or the anterior wall of the myocardium in standard parasternal short and long axis views. However, to apply 3D-BTI to other tissues in which structures can be oriented along the z-axis in sagittal and coronal views, tilted synthetic foci would be required. This could be achieved using subapertures in reception [START_REF] Tanter | Ultrafast compound imaging for 2-D motion vector estimation: application to transient elastography[END_REF] , which is the object of on-going work in our group.

Moreover, while the ultrasound acquisitions required to perform 3D-BTI are very short (~10 ms), postprocessing operations required up to twenty minutes to achieve the vectorial representation of a full volume. Further optimization of the algorithms along with the rapid growth of computational power of video cards and processors could potentially lead to a real-time implementation of 3D-BTI in the future. Segmentation of the epicardium and endocardium remained also challenging using the Bmode images provided by the BTI acquisition and required to be manually performed by a trained cardiologist. A dedicated anatomical high quality imaging sequence based on harmonic imaging could improve considerably the detection of the wall boundaries.

BTI may be affected by low ultrasound signal to noise ratio at large depths. The performance of BTI need to be investigated in more details for low SNR configuration. We can anticipate BTI to be robust to incoherent noise such as clutter noise and thermal noise, thanks to the suppression of a large part of incoherent noises by the coherence function. However, we can also anticipate that the performance of BTI will drop considerably in very low SNR situations where incoherent noise becomes higher that the coherent signals.

Finally, while presenting similar advantages in terms of portability, safety, and real-time imaging capabilities, the 3D ultrafast ultrasound system that we used is a unique laboratory prototype that is not currently available in clinical practice. In summary, 3D-BTI may constitute a unique tool for non-invasively evaluating the myocardial fiber structure of a patient suffering from myocardial fibrosis or hypertrophic cardiomyopathy for diagnosis, treatment monitoring and follow-up.

Materials and methods

Spatial coherence on a matrix array transducer

Ultrasonic spatial coherence was assessed experimentally by focusing an ultrasound wave in a medium and by receiving the backscattered echoes on a 2D matrix array transducer (Figure 7.A and B). The 2D spatial coherence function 𝑅(∆𝑥, ∆𝑦) was obtained by computing the auto-correlation of the signals received by pairs of elements i and j distant by ∆𝑥 and ∆𝑦, along the main two coordinate axis of the matrix array: Where 𝑆 𝑘 (𝑥 𝑘 , 𝑦 𝑘 , 𝑡) is the signal received on the element k of the matrix transducer with coordinates 𝑥 𝑘 , 𝑦 𝑘 after applying a time-delay to compensate the propagation path. ∆𝑥 and ∆𝑦 are the distances between and element i and element j:

𝑅(∆𝑥, ∆𝑦) = 1 𝑁 𝑥 -|∆𝑥|
∆𝑥 = 𝑥 𝑖 -𝑥 𝑗 ; ∆𝑦 = 𝑦 𝑖 -𝑦 𝑗 (2)
T1 and T2 represents the averaging temporal window (equivalent to depth).

In random media, i.e. in absence of fibers, ultrasound scatterers are randomly distributed, (figure 7.A), so that the coherence function is predicted by the so-called Van Cittert-Zernike elements that composed the pair. C represents the theoretical coherence function for a focused emission in a random media that we would obtain with a matrix probe of 𝑁 𝑥 × 𝑁 𝑦 elements.

Data acquisition and signal processing

In order to focus the ultrasonic wave at specific locations in the medium, conventional ultrasound imaging relies on applying time delays to the transmitted wave (figure 7). A large number of transmitted waves focused at different locations of the medium is then required to obtain one single image of the medium. In volumetric imaging, the number of transmitted waves can reach several thousands of waves which results in very low volume rates. To overcome this issue, our approach consisted in using ultrafast plane wave imaging with a coherent compounding approach [START_REF] Provost | 3D ultrafast ultrasound imaging in vivo[END_REF] . Coherent compounding can generate, at any location of the region of interest, synthetic focal zones using only a few tens of emissions. Therefore, with this approach, large volumes can be imaged at very high volume rate between 50 and 5000 volumes/s.

The emissions consisted in the transmission of several tilted 2D plane waves (i.e. figure 11.A).

Each plane wave was defined by two angles. All pairs of angles ranging from -6° to 6° with a step angle of 2° were emitted, for a total of 49 tilted plane waves. Plane waves were emitted from a 2D matrix array probe with 32 (x-axis) by 35 (y-axis) elements (3MHz, 0.3 mm pitch, 50% bandwidth at -3dB, Vermon) connected to a customized, programmable, 1024 channel ultrasound system described in Provost et al. [START_REF] Provost | 3D ultrafast ultrasound imaging in vivo[END_REF] . This system was designed and built in-house for 3D ultrafast imaging applications. The 1024 independent channels could be used simultaneously in transmission, whereas receive channels were multiplexed to 1 of 2 transducer elements. Therefore, each emission was repeated twice, with the first half of the elements receiving during the first emission, and the second half of the elements receiving during the second emission.

The received signals were sampled at 12 MHz, for all the 1024 elements and the radiofrequency data were recorded. The 3D images were computed off-line using a dynamic receive The coherence functions were averaged at each depth of the 3D volume using a lateral sliding window with a kernel of (0.15x0.15 mm).

An elliptic fit was applied on the averaged coherence functions on arbitrarily fixed 0.4 autocorrelation coefficient threshold. The direction of the main axis of the ellipse was used to define the direction of the fibers at each point of space (i.e. figure 11.D). Finally, a vector representation was used to visualize the fibers directions in 3D using the Amira software (Visualization Sciences Group, Burlington, MA) (i.e. figure 11.E).

Experimental Setup

1) In vitro experiments and histology

In vitro experiments were performed on five porcine fresh left ventricle myocardial samples embedded in a 2% gelatin gel. A cross-marker whose orientation matched the local cardiac circumferential-longitudinal coordinates was labeled on the myocardial sample. The probe was positioned at 15 mm of the epicardium. Each myocardial region examined by BTI was dissected from the intact porcine heart into a rectangular block (10x10x25 mm). This block was fixed in formalin for 48h. To keep the orientation of the tissue block in the wright direction, the basal aspect of the tissue block was labeled with tattoo ink. Then the tissue block was cut in 4 to 5 tissue slices which were embedded in paraffin. Each paraffin block was completely cut at 5 µm of thickness in serial sections. One histological section out of 50 was stained with H&E: the resulting gap between two stained adjacent sections was 250 µm. All the sections were scanned with a Hamamatsu Nano Zoomer (Hamamatsu City, Japan). Fiber angles in all the histological digitalized images studied were computed using the Hough transform [START_REF] Duda | Use of the Hough Transformation to Detect Lines and Curves in Pictures[END_REF] .

2 elements that composed the pair. C represents the theoretical coherence function for a focused emission in a random media that we would obtain with a matrix probe of 𝑁 𝑥 × 𝑁 𝑦 elements.

  figure 1 3D-BTI principle. A. Tilted plane waves are emitted and the associated

Figure 3 . 2 3D-

 32 figure 2

figure 4 5 Fiber

 45 figure 4

theorem 31 : 7

 317 figure 7

  focusing beamforming algorithm followed by coherent compounding (see figure11.B). The 3D volume had a lateral size of 9.6x9.6 mm and a 3 cm to 6 cm depth depending on the region of interest. The lateral sampling of the image was 0.3 mm laterally (32x32 lines) and the axial sampling was 0.05 mm. 2D coherence functions were computed at each point of the 3D volume (figure 11.C) using equation(1). The temporal window [𝑻 𝟏 𝑻 𝟐 ], was set to 1.6 µs (i.e. five periods at 3 MHz).

)

  In vivo experiments Animal experiments were performed in an open chest sheep on the mid-anterior region. The experimental procedure was approved prior to use by the Institutional Animal Care and Use committee of Hôpital Européen Georges Pompidou (PARCC) according to the European Commission guiding principles (2010/63/EU). 3D-BTI was thus performed within an entire cardiac cycle. The 3D-BTI acquisition was performed at a volume rate of 90 volumes/s during 1.2 seconds in order to cover more than one cardiac cycle. It enabled the reconstruction of 100 fiber volumes. During the experiment, the electrocardiogram was recorded. Human experiments were performed by a trained cardiologist. The volunteer was positioned on the left lateral decubitus position. The probe was placed in the parasternal view and 2D realtime imaging was performed to position the probe and image the antero-septal wall. The 3D BTI acquisition was then launch to acquire 1s of data at a rate 90 volumes/s. The study received the approval of the CPP (Paris Ile de France VI, N° 15-15).

figure 1 3D- 4 5

 145 figure 1 3D-BTI principle. A. Tilted plane waves are emitted and the associated

figure 6 7

 67 figure 6Transthoracic imaging of myocardial fiber orientation in the human heart.

  𝑥 𝑖 , 𝑦 𝑖 , 𝑡) 𝑆 𝑗 (𝑥 𝑗 , 𝑦 𝑗 , 𝑡)

	1 𝑁 𝑦 -|∆𝑦|	𝑗 ∑ ∑ 𝑖	√∑ ∑ 𝑇 2 𝑡=𝑇 1 𝑆 𝑖 (𝑇 2 𝑆 𝑖 (𝑥 𝑖 , 𝑦 𝑖 , 𝑡) 2 𝑆 𝑗 (𝑥 𝑗 , 𝑦 𝑗 , 𝑡) 2 𝑡=𝑇 1	(1)

Acknowledgements

This work was supported by the European Research Council under the European Union's Seventh Framework Program (FP/2007-2013) / ERC Grant Agreement n°311025 and by LABEX WIFI (Laboratory of Excellence ANR-10-LABX-24) within the French Program "Investments for the Future" under reference ANR-10-IDEX-0001-02 PSL. J.P. was funded in part by a Marie Curie International Incoming Fellowship under the Seventh Framework Program for research (FP7) of the European Union.

Author contribution statement

C.P, J.P, O.V and M.P wrote the main manuscript. C.P, V.F, J.P, J.L.G, O.V and M.P performed experiments. P.B. performed histology. M.T, M.F and M.P designed study. All authors reviewed the manuscript.

Additional Information

Competing financial interests : The authors declare no competing financial interests.


Sandra Dérozier

Robert Bossy

Louise Deléger

Mouhamadou Ba

Estelle Chaix

Valentin Loux

Hélène Falentin

Claire Nédellec

Omnicrobe, an open-access database of microbial habitats, phenotypes and uses extracted from text

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Omnicrobe, an open-access database of microbial habitats, phenotypes and uses extracted from text

The drastic increase in microbe descriptions, habitats, phenotypes and uses in databases, reports and papers presents a two-fold challenge for the access to the information: (1) a standard representation to integrate heterogeneous data and ( 2 c) normalization with taxa from the NCBI taxonomy [START_REF] Schoch | NCBI Taxonomy: a comprehensive update on curation, resources and tools[END_REF] and concepts from the OntoBiotope ontology [START_REF] Chaix | Text mining tools for extracting information about microbial biodiversity in food[END_REF].

d) relation extraction between entities.

All information is aggregated, merged and indexed to be accessible from the information system. Omnicrobe offers powerful ways to express simple and complex ontology-based queries.

Omnicrobe also exposes an API (Application Programming Interface) that allows users to automatically integrate microbe biodiversity knowledge in external information systems. Omnicrobe has been used to quickly target useful strains in a food innovation application [START_REF] Harlé | Diversity of the metabolic profiles of a broad range of lactic acid bacteria in soy juice fermentation[END_REF]. The Omnicrobe database was used to identify species capable of fermenting soybean juice. Out of 10 species represented by 229 strains tested in vitro, 9 species, i.e. 179 strains, acidified (pH<6) soybean juice within 48h.

In future work, we plan to extend Omnicrobe to include new data sources such as GOLD, GBIF, …

  ) the normalization of textual descriptions by semantic analysis. Recent information extraction technologies from the text mining domain offer powerful ways to detect and structure textual information along ontology-based representations. The Omnicrobe database contains around 1 million descriptions of microbe properties and is populated by an Information Extraction workflow: a) 6 information sources: -biological resource catalogues: Inrae CIRM, BacDive by DSMZ, -sequence database: GenBank, -scientific literature: PubMed abstracts in microbiology. b) entity recognition: microbe taxa, habitats, phenotypes and uses.

Figure 1 .

 1 Figure 1. Text mining process.

Figure 2 .

 2 Figure 2. Omnicrobe web interface.

Université Paris-Saclay, INRAE, MaIAGE, Jouy-en-Josas, France

Université Paris-Saclay, INRAE, BioinfOmics, MIGALE Bioinformatics Facility, Jouy-en-Josas 78350, France

INRAE, STLO, Rennes, France