Chuanyu Zhang

Xiaofeng Guo
email: xiaofeng.guo@esiee.fr

Philippe Brunet
email: philippe.brunet@univ-paris-diderot.fr

† Maxime Costalonga

Laurent Royon

Acoustic streaming near a sharp structure and its mixing performance characterization

Keywords: Acoustic Streaming, Sharp structure, Mixing performance, Fluorescence flow visualization

Acoustic streaming can be generated in microchannels by low frequency acoustic transducer in the vicinity of sharp structures. Close to the tip, the strong curvature induces bent trajectories on the time-periodic acoustic flow, locally enhancing the streaming-generating force. In this study, we investigate the influence of the sharp structure and vibration velocity on the streaming flow. The vibration velocities are characterized by directly visualizing the displacement of tracing particles and the generated acoustic streaming is observed using particle image velocimetry, under various operating conditions. By measuring the concentration of a fluorescence dye, we evaluate the mixing performance for different values of tip angle, vibration amplitude and flow rate through the microchannel. Our results confirm that intense streaming is generated under low-frequency (2.5 kHz) acoustic condition when the local curvature of the boundary is close to or smaller than the viscous boundary layer thickness. It is shown that the sharpest the edge tip, the largest the vortices

This Author Accepted Manuscript is a PDF file of an unedited peer-reviewed manuscript that has been accepted for publication but has not been copyedited or corrected. The official version of record that is published in the journal is kept up to date and so may therefore differ from this version.

Terms of use and reuse: academic research for non-commercial purposes, see here for full terms. https://www.springer.com/aam-terms-v1

Introduction

Increasing interests in milli-or micro-fluidic devices represent both advantages and challenges. Their advantages in efficiency and compactness make them attractive for a large range of industrial [START_REF] Whitesides | The origins and the future of microfluidics[END_REF], biological and pharmaceutical [START_REF] Elvira | The past, present and potential for microfluidic reactor technology in chemical synthesis[END_REF] applications. However, due to geometry and pressure drop constraints, the flow regime in most microfluidic devices is laminar. Consequently, transport phenomena are generally considered as a challenge. Acoustic Streaming (AS) has gained interest in recent years due to its ability to generate flows even at very low Reynolds number, therefore making it advantageous in various applications [START_REF] Gopinath | Convective heat transfer due to acoustic streaming across the ends of a Kundt tube[END_REF][START_REF] Boluriaan | Acoustic streaming: from Rayleigh to today[END_REF][START_REF] Whitesides | The origins and the future of microfluidics[END_REF][START_REF] Lee | Smart Cooling Technology Utilizing Acoustic Streaming[END_REF][START_REF] Legay | Improvement of heat transfer by means of ultrasound: Application to a double-tube heat exchanger[END_REF]Huang et al. 2013bHuang et al. , 2018b;;[START_REF] Elvira | The past, present and potential for microfluidic reactor technology in chemical synthesis[END_REF].

Physically, the AS phenomenon is generated in a quiescent fluid subjected to a periodical acoustic or mechanical vibration field. A net mean flow originates from the dissipation of the

A u t h o r a c c e p t e d m a n u s c r i p t

vibration energy by viscosity of the fluid, the latter can be either Newtonian or non-Newtonian. Theoretically, the phenomenon can be modelled by the time-averaged nonlinear term of the Navier-Stokes (NS) equation coupled with vibration perturbation. The first AS phenomenon was observed as early as 150 years ago but continuous research efforts have been done until today due to the diversity of mechanisms and applications [START_REF] Faraday | On a Peculiar Class of Acoustical Figures; and on Certain Forms Assumed by Groups of Particles upon Vibrating Elastic Surfaces[END_REF][START_REF] Eckart | Vortices and Streams Caused by Sound Waves[END_REF][START_REF] Nyborg | Acoustic Streaming due to Attenuated Plane Waves[END_REF][START_REF] Stuart | Double boundary layers in oscillatory viscous flow[END_REF][START_REF] Lighthill | ACOUSTIC STREAMINGt[END_REF][START_REF] Boluriaan | Acoustic streaming: from Rayleigh to today[END_REF].

Recently, the observation of intense flow and vortices generated by acoustic fields within sawtooth-shaped microchannels have brought a renewal of interest to boundary-driven streaming flows in microfluidics (Huang et al. 2013(Huang et al. , 2018;;Nama et al. 2014;[START_REF] Ovchinnikov | Acoustic streaming of a sharp edge[END_REF][START_REF] Tobias | Forces acting on a small particle in an acoustical field in a thermoviscous fluid Forces acting on a small particle in an acoustical field in a thermoviscous fluid[END_REF]. (Huang et al. 2013a) firstly try to apply such type of streaming into mixing enhancement. [START_REF] Ovchinnikov | Acoustic streaming of a sharp edge[END_REF]) conducted theoretical and experimental works to fundamentally investigate such a phenomenon. In their validation experiments, the sharp edge (tip of a needle) provides a 461 Hz mechanical vibration, i.e. the actuation inside the liquid. They pointed out that such streaming velocities can be of several orders of magnitude higher than Rayleigh streaming velocities under similar amplitudes of vibration.

They attributed the origin of this new type of streaming to the centrifugal force around the sharp vibrating structure. For instance, in eye surgery process, the jets coming from ultrasonically vibrating needles can reach a speed of as high as 2 m/s (Zacharias and Ohl 2013b). In addition to high streaming velocity, another significant characteristic of such streaming is its relatively low acoustic frequency (Huang et al. 2013a;[START_REF] Ovchinnikov | Acoustic streaming of a sharp edge[END_REF] and energy input, which can be particularly adapted to transfer enhancement in microfluidics geometries. In the realm of acoustic excitations on microchannels, the vast majority of solutions [START_REF] Bruus | Acoustofluidics 2: Perturbation theory and ultrasound resonance modes[END_REF]) uses waves of a few MHz or higher frequency, so that the acoustic wavelength in water c = c / f, with c=1430 m/s), is of the order of the largest

A u t h o r a c c e p t e d m a n u s c r i p t

dimension of the channel, namely around half a millimetre or so. The most common setups use two transducers facing each other or a transducer and a reflector [START_REF] Wiklund | Acoustofluidics 14: Applications of acoustic streaming in microfluidic devices[END_REF][START_REF] Barnkob | Acoustic radiation-and streaminginduced microparticle velocities determined by microparticle image velocimetry in an ultrasound symmetry plane[END_REF][START_REF] Lei | Acoustic streaming in the transducer plane in ultrasonic particle manipulation devices[END_REF], in order to generate a resonance condition when the channel width is equal to the half of c . The counterpart of this technique is the necessity to use expensive devices like high frequency transducers or amplifiers. On the contrary, the use of sharp structures is adapted to low-frequency acoustic waves, even below the ultrasound range, typically of a few kHz, allowing the use of low-cost transducers and amplifiers. While the acoustic wavelength, of several tens of cm, is then much larger than the channel dimensions, the sharp structures induce a local intensification of the acoustic actuation to the fluid.

With the advantages of high streaming velocity and low excitation frequency, the sharp-edge acoustic streaming could potentially efficiently improve processes such as mixing, particles selection, heat transfer, etc. Researches, especially those from Huang's group, reported the characterization of sharp-edge streaming flow and its application in mixing (Huang et al. 2013b(Huang et al. , 2018b;;Nama et al. 2014a), bio-particle and liquid-drop control [START_REF] Leibacher | Acoustophoretic cell and particle trapping on microfluidic sharp edges[END_REF][START_REF] Cao | A Microfluidic Device with Integrated Sonication and Immunoprecipitation for Sensitive Epigenetic Assays[END_REF][START_REF] Ozcelik | Acoustofluidic Rotational Manipulation of Cells and Organisms Using Oscillating Solid Structures[END_REF] as well as various on-chip devices [START_REF] Huang | A reliable and programmable acoustofluidic pump powered by oscillating sharp-edge structures[END_REF](Huang et al. , 2015b, a;, a;[START_REF] Bachman | Acoustofluidic devices controlled by cell phones[END_REF]Zhao et al. 2019).

Although experiments and simulations on acoustic streaming in sawtooth-shaped microchannel have been done in the aforementioned studies (Huang et al. 2013b;[START_REF] Ovchinnikov | Acoustic streaming of a sharp edge[END_REF], several issues remain to be investigated. In particular, the vibration conditions required to generate significant streaming need more quantitative measurements. The available literature only provides the electrical input into transducers (mainly excitation frequency and input voltage) which isn't a direct physical parameter to understand the generation of streaming. Other experimental variables such as clamping,

A u t h o r a c c e p t e d m a n u s c r i p t

transducer type and quality can make reproducing challenging. As a result, the physical interpretation of the fluid-wall interaction in a confined mini/micro-fluidic structure under acoustic excitation, as well as its influence on acoustic streaming generation and mass transfer performances is of both fundamental scientific interest and applicational significance.

In this study, we investigate the streaming flow properties and mixing efficiency of acoustic streaming generated in Y-shaped microchannels with sharp edged structures. We design an experimental setup that allows to visualize streaming vortices by high-speed imaging and to investigate the key influencing parameters such as acoustic amplitudes, sharp edge angle, as well as flow rate. A visualization based on fluorescence dye concentration tracking is used to assess the mixing performance of the Y-mixer. The paper is organized as follows: section 2 summarizes the AS theory near sharp edges, section 3 presents the experimental setup and the typical flow visualization; then, sections 4 and 5 present experimental results respectively on the flow properties and mixing performances, and finally section 6 concludes on the main findings as well as some perspectives.

Theoretical background

The usual framework of AS studies describes the velocity and pressure fields via a decomposition into a first-order time-periodic and a second-order steady components, with the acoustic vibration velocity (and the subsequent induced periodic flow) being of much larger magnitude than the steady (streaming) flow response [START_REF] Lighthill | ACOUSTIC STREAMINGt[END_REF][START_REF] Rayleigh | The Theory of Sound, Volume One[END_REF]).

In the time-dependent acoustic vibration field, velocity and pressure are denoted as:

A u t h o r a c c e p t e d m a n u s c r i p t

Eq. 1 where, , are the velocity and pressure variation amplitudes due to the acoustic excitation; , are respectively the first-order acoustic velocity and pressure at time t. Re means the real part of the complex terms.

The total velocity and pressure fields comprise thus two components:

Eq. 2
Where, , are the total velocity and pressure of fluid; , are the steady flow velocity and pressure.

The sharp-edge AS in this study uses acoustic wavelength ( c = c / f, in the order of ~m) much higher than the characteristic flow dimensions (the height of the sharp edge or width of the channel, <1 mm). Both acoustic and streaming velocities are much lower than the sound speed (c=1430 m/s in water), the flow in the present study can be treated as incompressible flow. Inserting the Eq. 2 into the incompressible Navier-Stokes and continuity equations, leads to two relationships with time-dependent terms in Eq. 3 and Eq. 4:

Eq. 3

Eq. 4

Averaging these two equations over a period yields to:

A u t h o r a c c e p t e d m a n u s c r i p t

Eq. 5 Eq. 6

Here, is the time averaged inertia term due to the first-order acoustic field, and is treated as a steady force term in the second-order steady equation.

Applying a spatial derivative to the first-order acoustic equations Eq. 3 and Eq. 4, we further reduce them to their equivalent form Eq. 7, which is the classical wave equation for the acoustic field:

Eq. 7

From the form of , it is clear that a homogeneous acoustic field cannot generate any streaming flow, as the effective force would be null. However, the sharp tip structures with strong local curvature on the channel walls induce strong heterogeneity in the acoustic oscillating flow in the fluid. As sketched in Figure 1, at a given time, the vibration of the fluid is uniformly distributed except for the local zone around the sharp edge, the latter being the source of AS. More specifically, close to the tip, both the orientation of the acoustic field and the vibration amplitude provide favorable condition of an intense F s very near the tip. Far from the tip, the force is null or negligible. Therefore, the non-uniformity of the acoustic field caused by the sharp structure makes acoustic streaming at relatively low frequencies (several kHz) possible. This is furtherly confirmed by our direct observation presented later in Figure 5 as well as in the Supplementary video (Supplementary material 1).

With acoustic frequencies of the order of kHz and low-viscosity liquids like water, the dissipation length is much larger than the system size, and hence most of the dissipation is due to the viscous friction at the solid-liquid interfaces. This friction is located in a viscous boundary layer of width . [START_REF] Ovchinnikov | Acoustic streaming of a sharp edge[END_REF]) point out that the streaming only Since we can briefly describe that the fluid is jetted from the tip, the resulting vortices and the strong directional flow will actively bring disturbance into the flow along the channel.

A u t h o r a c c e p t e d m a n u s c r i p t

The ensemble of the phenomena is potentially applicable in a large number of engineering process, such as mixing, heat transfer enhancement, etc. The piezoelectric transducer (diameter 35 mm and thickness 0.51 mm) delivers acoustic vibrations to the glass slide and to the whole channel stuck onto it, at various resonance frequencies from about 1 kHz up to 40 kHz. We chose to operate at one of these resonance peaks f, namely that at f=2.5 kHz. It turns out that the best operating conditions in terms of streaming flow were obtained at this frequency.

The Y-shaped Polydimethylsiloxane (PDMS) channel is designed by 2D photo-lithography.

First, a mould made of SU8 negative photoresist was fabricated: a 50 microns-thick SU8 liquid layer was spin-coated on a silicon wafer. After a soft baking at 65°C, the resist was exposed to UV through a photomask. The resist was then immersed in a developer bath and then hard-baked at 95°C. Then, PDMS (Sylgard 184) was thoroughly mixed with 10% in mass of curing agent, degassed in a vacuum chamber, and poured on the SU8 mould in order to constitute a 2.5 mm thick layer of PDMS mixture on top of the wafer. The whole is put in a stove at 65°C during 4 hours.

A u t h o r a c c e p t e d m a n u s c r i p t

The PDMS is then sealed and stuck on a glass microscope slide after a 1 min oxygen plasma treatment of both sides. The microchannel of a height of 50 µm is then formed between the PDMS and the glass slide. The plasma treatment enables the microchannel to withstand the pressure from inlets without leakage. Key geometrical dimensions of the Y-mixer are detailed in Figure 3b. Sharp edges with different angles (30º, 60º, 80º and 90º) were fabricated from various moulds. Hence, all channels were sealed on a similar glass slide, and rapid cure epoxy resist was used to permanently glue the transducers on each glass slide. It turns out that the best acoustic coupling was obtained this way. The radius of curvature of the different tips were measured with a 120X microscope with five repetition measurements of each angle. According to our measurement of the curvature diameter of the sharp edge, its value ranges from 2.8 µm to 10.3 µm, which are shown in Figure 4 and Table 1.

These values are all smaller than the viscous boundary layer thickness 11. The transducers are excited with periodic sinusoidal signal. The available range of peak-topeak voltage is between 0 V and V (up to V). It is important to notice that the speed of sound ( ) in water is 1 430 m/s and that in glass is 4 540 m/s, so the wavelengths of vibration ( ) at 2.5 kHz are respectively 0.57 m and 1.82 m. These values are far greater than the scales of the channel, whose width is 0.5 mm and length is 25 mm. The amplitude of the acoustic wave is therefore expected to be uniform in the whole channel.

To visualize the flow, the fluid is seeded with fluorescent particles (green polystyrene microspheres, diameter 4.9 m, Thermo Scientific). However due to the limited sensitivity of the high-speed camera, the best contrast and image quality is obtained by direct lighting Since the streaming velocity near the tip can be of comparable magnitude to the vibration velocity, which is different from most previous studies (Squires and Quake 2005) and will be discussed later (Section 4.3.2), another challenge is to observe and measure these two velocities separately. The following methods are used: first, to detect the motion of the particles during one vibration period and compute the streaming velocity, the frame rate of the camera is set to as high as 25 kHz and the exposure time as short as 25 µs. Under this condition, any individual particle appears as a small sphere of a few pixels because the exposure time is much smaller than the vibration period ( . A total of ten pictures are captured for every vibrating cycle, making the visualization of the acoustic vibration possible. To obtain the streaming velocity, we choose sequential images with a time interval of one acoustic period to "hide" the vibrational velocity. Then, to measure the acoustic vibration velocity , the frame rate of the camera is reduced to 2.5 kHz with an exposure time extended to 400 µs. Under this condition, particles in these pictures appear as threads of fixed length that corresponds to the amplitude of its acoustic-driven motion this leads to an intense streaming force very near the tip, while this force is null or negligible far from the tip. Therefore, the non-uniformity of the acoustic field caused by the sharp structure makes acoustic streaming at relatively low frequencies (several kHz) possible.

Acoustic vibration speed

As the transmission from electrical transducer power to the solid glass then to the fluid is rather complex [START_REF] Uchida | Investigation of acoustic streaming excited by surface acoustic waves[END_REF], we here only experimentally characterize the resulted vibration. Therefore, the vibration amplitude and velocity induced by the transducer are measured directly on the glass slide and indirectly in the fluid, the former by a laser vibrometer and the second by particle visualization. Their comparison allows us to identify the relative oscillation velocity on the fluid-solid interface, where the acoustic streaming develops.

On the one hand, a laser vibrometer (OFV-505 model, Polytec) is used to directly measure the vibration on the glass slide -near the location of the channels and after the PDMS bonding. Within the uncertainties of measurements, the relationship between the prescribed voltage and the local vibration velocity , is found to be linear. We measured for instance = 2.2 mm/s, 4.4 mm/s and 6.9 mm/s corresponding respectively to input peak-to-peak voltages of 10, 20 and 30 V, at 2.5 kHz. In terms of acoustic amplitude , such that , and as , it corresponds to displacement of the order of 1 micron or smaller on the glass.

On the other hand, acoustic fluid displacement 2A is extracted from high-speed imaging, as already shown in Figure 5. Figure 6 shows the linear relation between induced vibration velocity of the particles inside the liquid and the peak-to-peak voltage of the piezo-

A u t h o r a c c e p t e d m a n u s c r i p t

transducer, which is consistent with [START_REF] Franke | The highly sensitive optical measurement of absolute SAW amplitudes for power flow analysis[END_REF]. The error of measuring the vibration velocity mainly comes from the pixel resolution and uncertainties of 2A is within 1 pixel (representing 0.75 μm) which corresponds to ±6.1 mm/s. For instance, our results show values of the particles acoustic vibration velocity =37.8 mm/s, 50.5 mm/s, 70.5 mm/s, 85 mm/s, and 101.7 mm/s at 2.5 kHz and corresponding to a peak-to-peak voltage of 10 V, 15 V, 20 V, 25 V and 30 V, which will be the main experimental conditions of the present study. For higher vibration amplitudes with transducer inputs between 30 V and 60 V, the vibration velocities are also shown in Figure 6. However, due to experimental limitations (mainly exposure time of the fast camera in case of very high streaming velocity), our experimental conditions are only limited to under 30 V.

Since only the relative fluid-wall vibration accounts for the acoustic streaming phenomenon, and considering the measured solid vibration is lower than the tenth of the liquid vibration, hereafter we use the particles visualization-based vibration speed as the reference, i.e., .

The following discussions are all based on this reference (Figure 6).

A u t h o r a c c e p t e d m a n u s c r i p t

Flow visualization

In this section, we keep the piezoelectric transducer input to 2.5 kHz and use three different vibration velocities: 37.8 mm/s, 70.5 mm/s and 101.7 mm/s. Different tip angles are used: 30º, 60º, 80º and 90º, so as to different channel flow rates = 0, 1, 4, 8, 12, 16 µL/min.

We focus our following discussions on how the acoustic streaming disturbs the main flow.

General characteristic of the flow

Acoustic perturbation to the fluid can be clearly observed from the velocity field given by the PIV. As shown in Figure 7a, under an acoustic vibration of 37.8 mm/s at 2.5 kHz, the fluid around the sharp edge is jetted from the tip through the transverse direction of the microchannel. Since no throughput is imposed in that case, two symmetrical vortices appear besides the sharp edge. Outside of the region near the sharp edge (shown in the red square in the upper-right corner of Figure 7a), however, no fluid disturbance is observed, confirming that AS originates from the geometrical singularity of the tip. As already discussed in Figure 2d, for the "round" structure with , there is no evident streaming appearing around the structure. This confirms the hypothesis given by the perturbation theory [START_REF] Ovchinnikov | Acoustic streaming of a sharp edge[END_REF], that only when the tip curvature is narrower than the acoustic boundary layer, that the acoustic streaming could appear.

Indeed considering , it is clear that the stress is stronger in areas where the spatial variations of are sharper. Hence, the streaming vortices around the sharp structure can be attributed to the periodical curved motion of the fluid around the tip. The

A u t h o r a c c e p t e d m a n u s c r i p t

streaming force generates a jet from the tip, perpendicular to the channel longitudinal direction, which in turns induces the aforementioned vortices pair.

A u t h o r a c c e p t e d m a n u s c r i p t

In Figure 7b, we trace the velocity distribution along the y-direction (following the thick white line shown in the upper right corner, y=0 is taken at the tip) based on the PIV results. It clearly shows a rapid increase of the AS velocity from zero to a maximum value in a distance as short as 0.025 mm. After the peak, the AS velocity reduces gradually as its momentum transfer to the main fluid in the channel. The sharp decrease of the streaming velocity is consistent with Ovchinnikov's estimation [START_REF] Ovchinnikov | Acoustic streaming of a sharp edge[END_REF]), the force is localized within a distance of less than from the tip; therefore the velocity along y-direction decreases locally and rapidly after reaching the maximum value. Under the above condition, the peak value (maximum AS velocity) is = 15.4 mm/s and it depends significantly on the excitation amplitude, sharp edge angle, etc. It should be noted that uncertainties of streaming velocity mainly come from two aspects: the first is the coupling of the acoustic vibration with the streaming, which primarily concerns the maximal velocity value ; and the second lies in the pixel resolution. For the first, the error of is proportional to the streaming velocity and is estimated to be ±3% of . For the second, the pixel resolution (1 pixel=0.75 μm) during two sequential images of a time slaps of 400 μs results in an error of ±0.9 mm/s.

A u t h o r a c c e p t e d m a n u s c r i p t AUTHOR ACCEPTED MANUSCRIPT

© 2019 Springer-Verlag GmbH Germany, part of Springer Nature. The red line shows qualitatively the velocity evolution profile which is in accordance with the simulation results by [START_REF] Ovchinnikov | Acoustic streaming of a sharp edge[END_REF].

21

Maximum acoustic induced streaming velocity

As a signature of the acoustic streaming intensity, we extract from each velocity field the transverse component of the maximum velocity along the tip axis as shown in Figure 7b.

The Figure 8 illustrates the influence of acoustic velocity and tip angles of the sharp edge on the maximum AS velocity . Sharper edges and stronger acoustic excitation lead to higher maximal velocity. Under our test conditions, the highest maximum AS velocity is 73 mm/s and it appears at = 101.7 mm/s (peak-to-peak voltage of 30 V) with the sharpest edge angle 30° (with 2r c = 2.8 microns). This confirms the two velocities are in the same order of magnitude, with a ratio of . Under lower acoustic excitation, the AS effect is lower and the ratio becomes ( = 37.8 mm/s, = 20 mm/s, for the 30° sharp edge). The fact that and are not linearly related is consistent with Ovchinnikov et al.'s theory [START_REF] Ovchinnikov | Acoustic streaming of a sharp edge[END_REF]. However, the framework of their work is , and hence, a quantitative agreement with our results should not be expected.

Furthermore, higher streaming velocities are found for structures with sharper tip and smaller angle. This is consistent with that sharper tips lead to stronger Rayleigh stress near the tip, for the same vibration velocity and confirms previously observed trends (Nama et al. 2014b).

A u t h o r a c c e p t e d m a n u s c r i p t

Let us note that the above conclusion is based on results . When , the streaming vortex is strongly influenced and a detailed visualization on the vortex shapes is necessary. The vortex size takes larger values for stronger vibration velocities and for a sharper tip angle (smaller diameter than curvature): this is consistent with the similar trends observed for the maximal streaming velocity. At =101.7 mm/s and =30º, the vortex size can reach 0.48 mm. As increases, the acoustic streaming spreads over the whole channel width ( approaches 1) for almost all tip angles tested. In such a situation, the vortex size is limited by the microchannel width. On the other hand, with low acoustic excitation, only sharp edge with smaller angle can help achieve big vortices: at =37.8 mm/s, the value of can reach 0.76 for the 30˚ tip but only 0.27 for the 80˚ one.

A u t h o r a c c e p t e d m a n u s c r i p t AUTHOR ACCEPTED MANUSCRIPT

© 2019 Springer-Verlag GmbH Germany, part of Springer Nature.

28

Figure 11 shows the variation of vortex size at both downstream and upstream of a 30˚ tip, under different flow rates. As the flow rate increases from null to 1 µL/min, the vortices size reduces drastically from 0.42 to 0.29 (upstream) and 0.26 (downstream). The vortices continue to shrink after 4 µL/min and they totally disappear when the flow rate is higher than 16 µL/min. In most cases, the upstream vortex is smaller than the downstream one, i.e., between 4 µL/min and 16 µL/min. However, results suggest that there is a threshold of flow rate around 2 µL/min under which the upstream vortex is even larger than the downstream one. A possible underlying reason is that the main stream velocity is lower than the upstream longitudinal AS velocity coming from the counter-flow direction. In such a situation, part of the fluid of the upstream vortex still flows back to the sharp edge, enhancing the vortex development. Otherwise, the fluid goes away with the main flow, thus reducing the vortices size. At higher throughput, comparatively, the downstream vortex varies relatively more slowly with flow rate and is subjected to less influence from main flow, presumably because the vortex is squeezed into the low-pressure area formed by the sharp structure. Although there still exists the upstream vortex at large throughput, it plays little role on the main flow.

From the point of view of transverse mass transfer in the channel, we introduce the disturbance extension to quantify the perturbation of the streaming on the main flow.

The identification of requires the definition of a critical streamline, as shown in the sub- in the active one. The is defined as the distance from the tip to the tangential limit of the critical streamline.

Relative value is also defined for the disturbance extension with respect to the transverse distance from the tip, noted as .

Eq. 8

where w=0.5 mm is the channel width, and =0.18 mm represents the tip height.

As shown in Figure 12, despite that the upstream vortex may disappear, the disturbance extension might not be zero. At high flow rates like 12 µL/min, the perturbation could reach 14% of the channel ( =0.14). To sum up, the flow perturbation is closely coupled with the flow rate and the disappearance of vortices does not mean no acoustic streaming occurs.

Channel throughput should be considered as a critical parameter that is coupled with Figure 13 shows how the normalized concentration (in Eq. 9, based on Beer-Lambert law)

across the transverse direction of the channel evolves at different locations. One location A1 before and two locations A2-A3 after the sharp edge are sampled. We observe that at A1 ahead of the sharp edge, the two fluids begin to mix with each other driven by the upstream vortex. At the sharp edge position, the mixing process is furtherly improved by the intense jet. After the sharp edge, the mixing process is no longer influenced by the acoustic streaming and therefore only mass diffusion contributes to mixing, as shown between A2 and A3 in Figure 13.

The enhancement of the mixing process around the sharp edge can be divided into two steps: i) the upstream vortex pushes some (high-concentration) fluid against the main flow, directly making it into contact with the other stream; ii) the fluid jetted by the sharp structure as well as the downstream vortex furtherly disrupt the mixing interfaces. To quantify the mixing, we introduce the mixing index along the flow direction ( Eq. 9 and Eq. 10).

;

Eq. 9

,

Eq. 10

where stands for the normalized concentration; is the gray value of pure water; is the gray value of the unmixed blue dye liquid; is the gray value of pixels. In Eq. 10, is the gray value of pixel; is the average gray value of N pixels. In this study, is the average gray value of pixels across the width of the channel at certain pixel position along flow direction.

A u t h o r a c c e p t e d m a n u s c r i p t

Errors of the Mixing Index in this study are mainly caused by dye deposition near the channel border both on the glass slide and on the PDMS. For every experimental condition, five repeated measures under the same experiment setup are conducted and the maximum error for Mixing Index is 0.02. As shown in Figure 14, without acoustic excitation and as a result of the sole diffusion, the mixing index slowly decreases. At =70.5 mm/s, the mixing index decreases from 0.8 to 0.6 while crossing the sharp structure and at =101.7 mm/s, it furtherly decreases from 0.78 to 0.4. Therefore, the stronger the acoustic excitation, the higher the mixing index obtained.

A u t h o r a c c e p t e d m a n u s c r i p t

It is obvious that the mixing process within the vicinity regions ahead and after the sharp structure are accelerated by the upstream and downstream vortices discussed above.

Actually, diffusive mixing happens before the sharp edge, so the Mixing index does not start from 1 in the figures.

It should be noted that the width of the channel can influence the calculation of the Mixing Index, so the values within the shadow (beside the sharp edge) in Figure 14 cannot be a reliable criteria of the mixing. Only values before and after sharp edge are directly comparable.

Figure 15 shows the mixing index as a function of different tip angles. Sharper angles clearly allow better mixing degree at the outlet since the mixing index drops from 0.7 to 0.52 for the 30° tip but only to 0.65 for the 80° one. This is consistent with the stronger streaming effect, in terms of vortex diameter and disturbance extent, with sharper tips. Conditions: , and .

Conditions: and

A u t h o r a c c e p t e d m a n u s c r i p t

The influence of the flow rate on the mixing performance is shown in Figure 16. When the inlet flow rate increases from 0.5 µL/min to 1.5 µL/min, the mixing index at the outlet turns from 0.38 to 0.6, thus reducing the mixing degree. In case of even higher flow rate at the same acoustic condition, the mixing continues to degrade. For example, with a flow rate of = 4 µL/min, the mixing index is only achieved to 0.68. Above 1.5 µL/min, the effect of acoustic streaming on the mixing performance can be considered negligible. The short fluid residence time in the channel [START_REF] Guo | Residence time distribution on flow characterisation of multichannel systems: Modelling and experimentation[END_REF] in cases of high flow rate seriously weakens the mixing efficiency at the outlet. To sum up, with low flow rate, one single acoustically excited sharp edge can be enough to accomplish a homogeneous mixing; however, under larger flow rates (particularly above Q s =1.5 µL/min for the current study), the disturbance brought by a single sharp edge is limited and several additional similar structures should be introduced along the microchannel to achieve good mixing performance.

Conclusions and perspectives

The main outcome of this study is to investigate the mechanism of AS by direct particle visualization and demonstrate the potentiality for mixing enhancement in microchannels via acoustic streaming around sharp structures. Phenomenologically similar as the Rayleigh streaming, but with strong local Reynolds stress enhancement near the tip, the sharp-edge acoustic streaming is viewed as a promising one since it enables strong streaming disturbance to the flow with low energy input and low-cost transducers and amplifiers.

Indeed, within flat channels, no streaming would be generated at such frequency and amplitude. We use particle visualization and PIV treatment to characterize the AS conditions of appearance, spatial structure, magnitude, and its disturbance level relative to the main shown to be inefficient for generating acoustic streaming and mixing.

Secondly, sharper edges and higher vibration amplitude give stronger streaming flows.

Under our test conditions, the highest maximum streaming velocity is 75 mm/s and it is obtained for =101.7 mm/s with the sharpest edge angle 30°. Meanwhile the corresponding vortex can almost cover the whole channel until reaching the opposite channel wall. In terms of shapes of vortices, we observe symmetrical counter-rotating vortices when flow rate is zero, and the upstream/downstream symmetry is broken as the flow rate is turned on. In the latter case, the upstream vortex disappears as is set beyond a relatively low value, while the downstream one remains within a larger range of . The perturbation distance can also be as large as the channel width under low flow rate and with the sharpest tip, even for moderate vibration magnitude. In further studies, the sharp edge structure could be redesigned in an asymmetrical way so that the downstream vortices could be restricted while improving the upstream one, the latter being favorable to the mixing efficiency. Another ongoing work is the evaluation of heat transfer [START_REF] Guo | Multi-channel heat exchanger-reactor using arborescent distributors: A characterization study of fluid distribution, heat exchange performance and exothermic reaction[END_REF], residence time distribution [START_REF] Guo | Residence time distribution on flow characterisation of multichannel systems: Modelling and experimentation[END_REF], and micromixing performance using chemical iodide-iodate reactions [START_REF] Guo | Mixing performance assessment of a multi-channel mini heat exchanger reactor with arborescent distributor and collector[END_REF][START_REF] Dong | Mixing and residence time distribution in ultrasonic microreactors[END_REF]. In particular, sharp edges structures are already available in a large number of thermo-fluidic devices, whose performance could potentially benefit from AS generation.

Examples include the micro pin fins for efficient heat sink design [START_REF] Chiu | The heat transfer characteristics of liquid cooling heat sink with micro pin fins[END_REF][START_REF] Ambreen | Effect of fin shape on the thermal performance of nanofluidcooled micro pin-fin heat sinks[END_REF], or ratchet-like microchannel Knudsen pumping thanks to thermal transpiration [START_REF] Chen | Analysis of flow induced by temperature fields in ratchet-like microchannels by Direct Simulation Monte Carlo[END_REF][START_REF] Wang | Numerical simulation of thermal edge flow in ratchet-like periodically patterned micro-channels[END_REF], etc. In summary, AS generated by acoustic fields near sharp structures would also be a heat transfer enhancement that uses both passive and active techniques.

  photolithography techniques.

Figure 1 :

 1 Figure 1: Sketch of acoustic oscillating flow around a tip edge

Figure 2 :Figure 3

 23 Figure 2: Mechanism of acoustic streaming under zero flowrate (a) acoustic streaming condition, sharp-edge curvature should be lower than the acoustic boundary layer 2 < ; (b) fluorescent particles distribution around a sharp-edge ( )without acoustic field; (c) the streamline of the acoustic excited particles movement under a stack of 100 images equivalent to a duration of 100 ms, = 37.8 mm/s, ; (d) round-edge with a curvature of 100 m under acoustic excitation but with no acoustic streaming, , = 37.8 mm/s, .

  5 µm, determined by (Boluriaan and Morris 2003), with the kinematic viscosity of the fluid, and the angular frequency.

Figure 3 :Figure 4 :

 34 Figure 3: Sketch of the experimental setup (a) fluid circulation and visualization system;(b) Y-type mini-channel with a sharp edge

  with white light, hence with the diffused light instead of the fluorescence-emitted light. To enhance the contrast, a cold-light beam shines from the bottom of the glass slide. The depth of field of the microscope lens is about 10 microns, five-times smaller than the channel depth (50 microns). Hence, after careful focal adjustments (based on finding the position that shows the maximal transverse velocity for a given flow rate), we obtain high fidelity pictures with clear focus, as shown in Figure2b-d.

(

  Figure 5). Results are shown in supplementary video. During the experiment, the pictures are continuously synchronized with the camera to the computer. The obtained images are then treated with open source software ImageJ (Schneider et al. 2012). For vibration velocity, the particles' displacements are directly measured from the A u t h o r a c c e p t e d m a n u s c r i p t captured images, and the results as well as the uncertainties are discussed in detail in Section 4.2. For the streaming field, the relative displacement of particles at a given phase of the vibration period and at given locations, allows the determination of the velocity field. Successive frames are converted into an array of displacement vectors for each particle by the software PIVlab (Thielicke and Stamhuis 2014). Second, summing the intensity of the images of a number of frames allows clear visualization of trajectories and then the location, form and sizes of vortices. Besides, COMSOL Multiphysics © simulation is used as a control study to verify the reliability of the PIV result in terms of velocity characterization. The control study is done without acoustic excitation using a standard laminar flow model. The results are shown to be highly close (See Supplementary material 2).

Figure 5 :

 5 Figure 5: The mechanism of acoustic streaming explained by high-speed visualization of acoustic vibrations within the fluid, from particle trajectories over a few periods, the velocity measured is 101.7 mm/s. (see Supplementary video for dynamic visualization effect).

Figure 6 :

 6 Figure 6: Acoustic vibration velocity measured by particles visualization at different input peak-to-peak voltage of the piezo-transducer, fitting line .

Figure 7 :

 7 Figure 7: Acoustic streaming velocity field, at and . (a) velocity field obtained from PIV and (b) velocity distribution along the transverse direction from the tip. The acoustic streaming is obtained from a tip with The red

Figure 8 :

 8 Figure 8: Maximum streaming velocity versus acoustic velocity, generated near sharp tips of different angles, under no-flow condition ( .

Figure 9 :Figure 10 :

 910 Figure 9: Flow visualization near the sharp tip, using particle trajectories at different channel throughputs, experimental conditions: , mm/s and . The flow direction is from left to right.

Figure 11 :

 11 Figure 11: Upstream and downstream vortex size (D v and D v /w) at different flowrates Condition: , and .

Figure 12 .

 12 Figure 12. The critical streamline (Wang et al. 2012) separates an active perturbation zone

Figure 12 :

 12 Figure 12: Disturbance extension D dis as a function of channel throughput Conditions: , and .

Figure 13 :

 13 Figure 13: Concentration field across the microchannel section in the nearby area before and after the sharp edge. X axis is the normalized concentration (value 0 and 1 represents two different fluids, 0.5 means complete mixing); Y axis is channel width. Conditions: , , and .

Figure 14 :Figure 15 :

 1415 Figure 14: Mixing Index at different acoustic excitation amplitudes. Conditions: , and .

Figure 16 :

 16 Figure 16: Mixing index at different flow rates, the values shown being those for each individual inlet.

  flow characteristics under different experimental conditions (tip angle, vibration amplitude/velocity, flow rate) are quantitatively investigated. Our mains conclusions are the following: Firstly, regarding the AS phenomenon, our study confirms the perturbation theory about the boundary layer condition. Sharp edges, whose curvature radius (1.3 to 5.2 m in this work) is smaller than the viscous boundary layer (11.3 m in this work), can achieve strong acoustic streaming flow. Conversely, a round edge structure with a curvature diameter of 100 m, is

Finally, perhaps not

  implausibly, stronger streaming (with sharper tip and under larger vibration magnitude) helps achieve better mixing. Significant mixing improvement is witnessed with the fluorescent dye measurement: a mixing index of 0.38 (mostly mixed) when AS is generated, compared to 0.73 when acoustic excitation is off, for the case of 2 µL/min, 30° tip and under a vibration amplitude of 70.5 mm/s. In terms of flowrate, mixing enhancement by streaming is obviously weakened at flowrates above 1.5 µL/min. In very intuitively, the disturbance extent is crucial for the mixing application, but not the vortex size itself.

  -Verlag GmbH Germany, part of Springer Nature. 43 Zacharias J, Ohl CD (2013b) Fluid dynamics, cavitation, and tip-to-tissue interaction of longitudinal and torsional ultrasound modes during phacoemulsification. J Cataract Refract Surg 39:611-616. doi: 10.1016/j.jcrs.2012.10.050 Zhao S, He W, Ma Z, et al (2019) On-chip stool liquefaction via acoustofluidics. Lab Chip 19:941-947. doi: 10.1039/C8LC01310A

© 2019 Springer-Verlag GmbH Germany, part of Springer Nature.

Acknowledgement

The authors would like to acknowledge the scholarship support from the CSC (China Scholarship Council). Eric Falcon is acknowledged for his help in the measurement of the vibration velocity with the laser vibrometer. They would like to thank the MSC and LIED lab members for their daily support and discussion.


Senior Member, IEEE Esa Ollila

Fellow, IEEE Daniel P Palomar

Senior Member, IEEE Frédéric Pascal

Affine equivariant Tyler's M-estimator applied to tail parameter learning of elliptical distributions

Keywords: Tyler's M-estimator, scatter matrix, covariance matrix, elliptical distributions

We propose estimating the scale parameter (mean of the eigenvalues) of the scatter matrix of an unspecified elliptically symmetric distribution using weights obtained by solving Tyler's M-estimator of the scatter matrix. The proposed Tyler's weightsbased estimate (TWE) of scale is then used to construct an affine equivariant Tyler's M-estimator as a weighted sample covariance matrix using normalized Tyler's weights. We then develop a unified framework for estimating the unknown tail parameter of the elliptical distribution (such as the degrees of freedom (d.o.f.) ν of the multivariate t (MVT) distribution). Using the proposed TWE of scale, a new robust estimate of the d.o.f. parameter of MVT distribution is proposed with excellent performance in heavy-tailed scenarios, outperforming other competing methods. R-package is available that implements the proposed method.

I. INTRODUCTION

W E MODEL the observed p-variate observations x 1 , . . . , x n as independent and identically distributed (i.i.d.) random samples from an unspecified centered (i.e., symmetric around the origin) elliptically symmetric (ES) distribution [START_REF] Fang | Symmetric Multivariate and Related Distributions[END_REF], [START_REF] Ollila | Complex elliptically symmetric distributions: survey, new results and applications[END_REF]. A continuous random vector x ∈ R p has centered ES distribution if it possesses a probability density function (pdf) of the form f (x) = C p,g |Σ| -1/2 g(x ⊤ Σ -1 x), where g : R ≥0 → R >0 is called the density generator, Σ ≻ 0 is the positive definite symmetric matrix parameter, called the scatter matrix, and C p,g is a normalizing constant ensuring that f (x) integrates to 1. We let x ∼ E p (0, Σ, g) to denote this case. For example, the centered multivariate normal (MVN) distribution N p (0, Σ) is obtained when g(t) = exp(-t/2) while the multivariate t (MVT) distribution with ν > 0 degrees of freedom (d.o.f.) is obtained when

g(t) = (1 + t/ν) -(p+ν)/2 . ( 1
)
Parameter ν > 0 is a tail parameter of the density. For ν → ∞, the MVT distribution reduces to the MVN distribution, while ν = 1 corresponds to the multivariate Cauchy distribution. Also, many other subclasses of ES distributions are parametrized by a density generator indexed by an additional tail parameter that is unknown in practice and needs to be estimated. Learning this unknown parameter is also one of the goals of this paper.

Esa Ollila is with the Department of Information and Communications Engineering, Aalto University, FI-00076 Aalto, Finland (e-mail: esa.ollila@aalto.fi). Daniel P. Palomar is with the Hong Kong University of Science and Technology, Hong Kong (e-mail: palomar@ust.hk). Frédéric Pascal is with the CNRS Université Paris-Saclay, CentraleSupélec, 91190 Gifsur-Yvette, France (e-mail:frederic.pascal@centralesupelec.fr).

We are interested in estimating the scale parameter of the scatter matrix, defined as the mean of its eigenvalues,

η = tr(Σ) p = 1 p p i=1 λ i , (2)
where

λ i > 0 denotes the ith eigenvalue of Σ. Formally, η ≡ η(Σ) is a scale parameter if it verifies η(I) = 1 and η(aΣ) = aη(Σ) for all a > 0 [3].
Tyler's M-estimator [START_REF] Tyler | A distribution-free M-estimator of multivariate scatter[END_REF] is a popular robust M-estimator of the scatter matrix that has been extensively studied both in signal processing and statistics literature (e.g, [START_REF] Sun | Regularized Tyler's scatter estimator: Existence, uniqueness, and algorithms[END_REF]- [START_REF] Romanov | Tyler's and Maronna's Mestimators: Non-asymptotic concentration results[END_REF]). Tyler's M-estimator is defined as the solution to the fixedpoint equation

Σ = 1 n n i=1 p x ⊤ i Σ-1 x i x i x ⊤ i H( Σ; {x i }). (3)
Note that me way also write the map H(•; •) in the form

H( Σ; {x i }) = 1 n n i=1 ŵi x i x ⊤ i with ŵi = p x ⊤ i Σ-1 x i ,
where ŵi , i = 1, . . . , n, are referred to as Tyler's weights.

Tyler's M-estimator is unique only up to a scaling factor, and therefore a common convention is to consider a solution that verifies tr( Σ) = p. Thus, Tyler's M-estimator is actually an estimator of a shape matrix (normalized scatter matrix) Λ, defined by Λ = Σ/η = pΣ/ tr(Σ), and verifying tr(Λ) = p.

In this paper, we propose an estimator of the scale η based on Tyler's weights ŵi . The proposed scale estimate along with Tyler's M-estimator Σ are then jointly used for constructing affine equivariant robust estimates of the scatter matrix Σ and the covariance matrix R = cov(x) (or their shrinkage versions). These developments are described in Section II. Then, in Section III, we propose a unified framework allowing to estimate the tail parameter of the elliptical distribution using the proposed scale statistic η. In the case of the MVT distribution, this leads to a new estimate of the d.o.f. parameter based on Tyler's weights. Finally, Section IV demonstrates the relevance of the proposed approach on simulated data, with concluding remarks in Section V. In the R package fitHeavyTail [START_REF] Palomar | fitHeavyTail: Mean and Covariance Matrix Estimation under Heavy Tails[END_REF], the function fit_Tyler implements this method.

II. ESTIMATE OF SCALE, SCATTER, AND COVARIANCE

MATRIX BASED ON TYLER'S WEIGHTS Assuming that x ∼ E p (0, Σ, g) has finite 2nd-order moments, then its covariance matrix,

R = E[xx ⊤ ] satisfies R = θ • Σ for θ = E[r 2 ] p , (4)
where r 2 = Σ -1/2 x 2 is the squared Mahalanobis distance of x w.r.t. Σ, whose pdf is given by

f r 2 (t) = Ct p/2-1 g(t), (5)
where C = ∞ 0 t p/2-1 g(t)dt. Hence pdf of r 2 has a oneto-one correspondence with density generator g. From (4), we notice that the scatter matrix Σ is proportional to the covariance matrix R (assuming R exists). In the MVN case, θ = 1, while for the MVT distribution with density generator as in (1) one obtains θ = ν/(ν -2) for all ν > 2.

A. Estimate of scale

As mentioned earlier, Tyler's M-estimator Σ loses information of the scale η. However, it is yet possible to construct an estimate of η from Tyler's weights. Our Tyler's weights-based estimate (TWE) of scale is defined as the harmonic mean of reciprocal of weights, 1/ ŵi 's, that is,

ηTWE = 1 n n i=1 ŵi -1 = p n n i=1 [x ⊤ i Σ-1 x i ] -1 -1 . (6)
One can esily verify that this is a proper scale estimate in the sense that if η * TWE is computed on scaled observations,

x * i = c • x i , i = 1, . . . , n, then η * TWE = c 2 • ηTWE . This follows because Tyler's M-estimator Σ in (3) with trace constraint is invariant to scaling the data, so Σ * = Σ.
The proposed estimate (6) can also be motivated from the following result derived in the high-dimensional random matrix theory (RMT) regime, where p, n → ∞ with n > p and their ratio tending to constant: p/n → c ∈ (0, 1). Namely, let Σ be Tyler's M-estimator in (3) verifying tr( Σ) = p. Then, it was shown in [START_REF] Zhang | Marčenko-pastur law for Tyler's M-estimator[END_REF], [START_REF] Romanov | Tyler's and Maronna's Mestimators: Non-asymptotic concentration results[END_REF] that max ℓ |η ŵℓ -1| → 0 almost surely. The authors in [START_REF] Zhang | Marčenko-pastur law for Tyler's M-estimator[END_REF] derived this result for the case that data is i.i.d. Gaussian N p (0, I) while [START_REF] Romanov | Tyler's and Maronna's Mestimators: Non-asymptotic concentration results[END_REF] extended these results for more general distributions. Thus since 1/ ŵi concentrate on η, (6) is a natural robust estimator of scale. Many other robust scale statistics could be constructed from 1/ ŵi 's, such as the median, trimmed mean, etc. In the next subsection, we illustrate why the proposed harmonic mean [START_REF] Pascal | Generalized robust shrinkage estimator and its application to STAP detection problem[END_REF] is the most natural.

B. Affine equivariant estimate of scatter matrix

Using the scale estimate ηTWE and Tyler's M-estimator Σ (with tr( Σ) = p), we can form an estimate of the scatter matrix

ΣTWE = ηTWE • Σ (7)
referred to as TWE of scatter matrix. Thus ηTWE is scale statistic derived from ΣTWE since ηTWE = tr( ΣTWE )/p. Equivalently, the trace of ΣTWE can be easily computed as the harmonic mean of Tyler's quadratic form:

tr( ΣTWE ) = 1 n n i=1 [x ⊤ i Σ-1 x i ] -1 -1
.

Recalling [START_REF] Pascal | Generalized robust shrinkage estimator and its application to STAP detection problem[END_REF] we can write [START_REF] Ollila | Regularized M -estimators of scatter matrix[END_REF] in the following more intuitive form:

ΣTWE = 1 n n i=1 vi x i x ⊤ i , (8a) vi = ŵi 1 n n ℓ=1 ŵℓ = [x ⊤ i Σ-1 x i ] -1 1 n n ℓ=1 [x ⊤ ℓ Σ-1 x ℓ ] -1 (8b)
where v1 , . . . , vn are normalized Tyler's weights that verify

1 n n i=1 vi = 1.
Eq. ( 8a) and (8b) illustrate that ΣTWE is a weighted sample covariance matrix (SCM) with weights vi .

Finally, we draw the parallel of Tyler's M-estimating equation and our estimator [START_REF] Chen | Robust shrinkage estimation of high-dimensional covariance matrices[END_REF]. First, note that Tyler's M-estimating equation ( 3) verifies broader invariance than just invariance with respect to scaling of the data matrix. Namely, denoting the unit-norm normalized data by xi = x i / x i , i = 1, . . . , n, one can easily verify that the fixed-point equation in (3) can be rewritten as Σ = H( Σ; {x i }), so based on normalized data. Furthermore, since tr( Σ) = p, one has:

tr( Σ) = 1 n n i=1 p x⊤ i Σ-1 xi tr xi x⊤ i = p, or equivalently 1 n n i=1 [x ⊤ i Σ-1 xi ] -1 = 1 since tr xi x⊤ i = 1.
It follows that Tyler's M-estimator Σ with tr( Σ) = p is the solution to the following fixed-point equation:

Σ = 1 n n i=1 wi xi x⊤ i with wi = p [x ⊤ i Σ-1 xi ] -1 1 n n ℓ=1 [x ⊤ ℓ Σ-1 xℓ ] -1 , (9)
where wℓ 's are the normalized Tyler's weights computed on normalized (unit norm) observations xi 's. Thus, while Tyler's M-estimator Σ can be interpreted as a weighted SCM based on normalized data {x i } as shown in [START_REF] Couillet | Large dimensional analysis and optimization of robust shrinkage covariance matrix estimators[END_REF], TWE of scatter ΣTWE in (8) can be viewed as weighted SCM of actual (nonnormalized) data {x i }.

It is worthwhile to point out that most robust estimators of scatter are affine equivariant in the sense that an affine transformation on the data x i → Ax i , i = 1, . . . , n, induces following transformation on the estimate:

Σ({Ax i }) = A Σ({x i })A ⊤ , ∀A ∈ R p×p invertible. ( 10
)
For example, robust Maronna's [START_REF] Maronna | Robust M-estimators of multivariate location and scatter[END_REF] M-estimators, Sestimators [START_REF] Davies | Asymptotic behaviour of S-estimates of multivariate location parameters and dispersion matrices[END_REF], or MM-estimators [START_REF] Tatsuoka | On the uniqueness of S-functionals and M-functionals under nonelliptical distributions[END_REF] are affine equivariant. However, Tyler's scatter matrix Σ is not affine equivariant since [START_REF] Soloveychik | Performance analysis of Tyler's covariance estimator[END_REF] only holds up to multiplicative scalar factor as shown in [START_REF] Zhang | Marčenko-pastur law for Tyler's M-estimator[END_REF]. Affine equivariance is desirable since if x ∼ E p (0, Σ, g), then Ax ∼ E p (0, AΣA ⊤ , g). Hence the transformed data shares the same elliptical distribution, but the scatter matrix parameter is mapped to AΣA ⊤ . Thus a natural requirement to be imposed on any scatter matrix estimator is that it should verify this same equivariance principle under transformations x i → Ax i . This is shown next.

Lemma 1. TWE of scatter matrix ΣTWE is affine equivariant: ΣTWE ({Ax

i }) = A Σ({x i })A ⊤ , ∀ invertible A ∈ R p×p .
Proof. It is straightforward to verify that Tyler's M-estimator (with tr( Σ) = p) is equivariant in the sense that if

Σ * = Σ({x * i }) denotes Tyler's M-estimator (verifying tr( Σ * ) = p) computed on data x * i = Ax i , i = 1, . . . , n, then Σ * = cA ΣA ⊤ , c = p tr(A ΣA ⊤ ) (11)
for all invertible A ∈ R p×p . Now let ŵ * i , i = 1, . . . , n denote the corresponding Tyler's weights. Then the scale estimate (6) computed as the harmonic mean of reciprocals of weights

1/ ŵ * i = [x * i ] ⊤ ( Σ * ) -1 x * i is η * TWE = 1 n n i=1 ŵ * i -1 = p n n i=1 [x ⊤ i A ⊤ ( Σ * ) -1 Ax i ] -1 -1 = c -1 p n n i=1 [x ⊤ i Σ-1 x i ] -1 -1 = c -1 ηTWE (12)
where in the 2nd to last identity we simply utilized [START_REF] Zhang | Marčenko-pastur law for Tyler's M-estimator[END_REF]. Thus it follows that

Σ * TWE = η * TWE • Σ * = c -1 ηTWE • cA ΣA ⊤ = A(η TWE Σ)A ⊤ = A ΣTWE A ⊤ .
In the case of large dimensional data, one can also consider a shrinkage Tyler's M-estimator of the scatter matrix as

ΣTWE,β = β ΣTWE + (1 -β)η TWE I, (13)
where the data adaptive shrinkage parameter β ∈ [0, 1] is computed as described in [18, Sect. IV.C]. However, unlike the estimator in [START_REF] Ollila | Shrinking the eigenvalues of M-estimators of covariance matrix[END_REF], the shrinkage TWE in [START_REF] Romanov | Tyler's and Maronna's Mestimators: Non-asymptotic concentration results[END_REF] provides an estimator of scatter instead of shape matrix.

C. An estimator of covariance matrix

If the density generator g (and hence the underlying ES distribution) is specified, then the value of θ in (4) can be determined, and we can use relationship (4) to obtain a covariance matrix estimator as RTWE = θ • ΣTWE . For example, if the data has an MVN distribution, then θ = 1 while θ = ν/(ν -2) in the case of an MVT distribution with ν d.o.f. However, often the underlying parametric family is known, but the underlying tail parameter, say ν, indexing the density generator is unknown. As is shown in Section III, we can form an estimate of ν, denoted νTWE , using Tyler's weights. Since θ = h(ν) (cf. Eq. ( 15) below), a TWE of covariance matrix can be computed as RTWE = θTWE • ΣTWE , where θTWE = h(ν TWE ).

III. ESTIMATING THE TAIL PARAMETER OF ES DISTRIBUTION

From (4) we can induce the following relationship between the scale parameter η cov = p -1 tr(R) of the covariance matrix and scale η = p -1 tr(Σ) of the scatter matrix:

η cov = θη ⇔ θ = η cov /η. ( 14
)
Note that a natural estimate of η cov is p -1 tr(S), where S = 1 n n i=1 x i x ⊤ i denotes the sample covariance matrix (SCM). On the other hand, if density generator g is specified up to unknown tail parameter ν, thus indexed by g ν (•), then θ in (4) is a following function of the tail parameter ν:

θ = ∞ 0 tf r 2 (t; ν)dt h(ν) (15)
where the pdf f r 2 (•), defined in [START_REF] Sun | Regularized Tyler's scatter estimator: Existence, uniqueness, and algorithms[END_REF], is one-to-one with g ν (•). We do not need numerical integration in most practical cases as often closed-form expression for h(ν) can be derived. Then, after solving the inverse mapping, ν = h -1 (θ) = h -1 (η cov /η), Algorithm 1 offers a unified approach for estimating the tail parameter of an ES distribution:

Algorithm 1: Distribution tail parameter learning Input:

Data {x i } n i=1
Output: Estimated tail parameter ν 1. Compute Tyler's M-estimator Σ and weights ŵi 's in (3); 2. Compute ηTWE in ( 6), and set θTWE = p -1 tr(S)/η TWE ; 3. Using [START_REF] Maronna | Robust M-estimators of multivariate location and scatter[END_REF], estimate ν as νTWE = h -1 ( θTWE ).

As an example, if x follows an MVT distribution with ν > 2 d.o.f., one has that θ = h(ν) = ν ν-2 , which unfolds the relation:

ν = h -1 (θ) = 2θ θ -1 for θ > 1. (16)
Note that ν > 2 is required for the covariance matrix R to exist. The obtained estimator νTWE is closely related to estimator in [18, Alg. 1], referred to as OPP estimator for short. OPP is an iterative approach that iteratively (re-)computes the maximum likelihood estimator (MLE) Σ of the MVT distribution with ν given by the current estimate of d.o.f. parameter ν (k) . It then computes θ = tr(S)/ tr( Σ) which provides an update ν (k+1) = h -1 ( θ) via [START_REF] Davies | Asymptotic behaviour of S-estimates of multivariate location parameters and dispersion matrices[END_REF]. The algorithm iterates for k = 0, 1, 2 . . . until convergence, starting from an initial start ν (0) = νkurt , where νkurt is an estimate of ν based on elliptical kurtosis, proposed in [START_REF] Ollila | Shrinking the eigenvalues of M-estimators of covariance matrix[END_REF], and referred to as kurtosis estimator. We also proposed an improved version of OPP estimator in [START_REF] Pascal | Improved estimation of the degree of freedom parameter of multivariate t-distribution[END_REF], which, however, is impractical for large n and p. ML estimation of ν via the Expectation-Maximization (EM) approach is considered in [START_REF] Liu | ML estimation of the t-distribution using EM and its extensions, ECM and ECME[END_REF]. This method is unfortunately rather unstable [START_REF] Fernandez | Multivariate Student-t regression models: Pitfalls and inference[END_REF].

IV. SIMULATION STUDIES

We first consider the case where scatter matrix Σ has an autoregressive model (AR(1)) structure, (Σ) ij = η̺ |i-j| , where η = tr(Σ)/p is the scale parameter and ̺ is the correlation parameter, ̺ ∈ (-1, 1). Since Tyler's M-estimator Σ is invariant to the data scaling, we can set η = 1 without favoring any estimator over the other. The number of Monte-Carlo runs is 5,000, and samples are generated from an MVT distribution with different choices of d.o.f. parameter ν.

First, we investigate how the TWE of d.o.f. parameter ν compares against OPP and kurtosis estimators. Figure 1 displays the boxplots in the case that p = 100 and varying sample lengths when ν = 5 or ν = 3. As can be noted, the proposed TWE attains the best accuracy as well as the smallest variability. Moreover, for ν = 5, its median values are right on the spot. The kurtosis estimator obviously performs poorly when ν = 3 since the 4th-order moment does not exist in this case. Figure 2 shows the average mean squared error (MSE), (ν -ν) 2 , which further illustrates the benefits and high accuracy of the proposed TWE against its competitors. Figure 3 displays the boxplots of different estimates of scale η. Here we compare ηTWE to OPP estimate of scale, defined as ηOPP = tr( Σ)/p, where Σ is the MLE of scatter based on ν = νOPP . We also compare with the scale estimate provided by the SCM, defined as ηcov = tr(S)/p, but multiplied by θ -1 = (ν -2)/ν to obtain an estimate of η; recall [START_REF] Palomar | fitHeavyTail: Mean and Covariance Matrix Estimation under Heavy Tails[END_REF]. We can notice from Figure 3 that TWE slightly underestimates the true scale while OPP is overestimating. We also notice that the SCM estimator is clearly unbiased for ν = 5, but has huge variability. Figure 4 displays the median values of ν for range of ν values when n = 150 (and p = 100 and ̺ = 0.6 as earlier). The proposed TWE estimator significantly outperforms the other estimators for all d.o.f. ν ∈ [START_REF] Paindaveine | A canonical definition of shape[END_REF][START_REF] Chen | Robust shrinkage estimation of high-dimensional covariance matrices[END_REF]. We now consider an example based on stock market data. We generate synthetic data (p = 100 assets) with heavy tails following MVT distribution with d.o.f. ν = 4 and covariance matrix as measured from stocks of S&P 500 index. Figure 5 compares the estimated value of ν versus the number of observations for the following methods: kurtosis estimator, OPP estimator [START_REF] Ollila | Shrinking the eigenvalues of M-estimators of covariance matrix[END_REF], and the proposed estimator, with the latter being clearly superior, illustrating its promising performance for real-world financial data.

V. CONCLUDING REMARKS

We proposed a new robust estimator of scale parameter of an elliptical distribution based on the weights from Tyler's M-estimator, which was further used to construct an affine equivariant Tyler's M-estimator. We then proposed a unified framework to estimate the tail parameter of an elliptical distribution. Finally, it should be noted that this method generalizes to complex-valued data in a straightforward manner.

Fig. 1 :Fig. 2 :Fig. 3 :

 123 Fig. 1: Boxplots d.o.f. estimates ν as a function of n when ν = 5 (top) and ν = 3 (bottom); p = 100, ρ = 0.6.

Fig. 4 :Fig. 5 :

 45 Fig. 4: Estimated values of ν; p = 100, ρ = 0.6, n = 150. The black line indicates the true value.