Francisco A R Pereira

Kaline S Sousa

Graycyelle R S Cavalcanti

Denise B França

Líbia N F Queiroga

Ieda M G Santos

Maria G Fonseca

Maguy Jaber

Green biosorbents based on chitosan-montmorillonite beads for anionic dye removal

Keywords: Clay minerals, biopolymer, composites, wastewater, organic pollutants

Highlights

 Chitosan-montmorillonite beads composites were synthesized  The key parameters that control the adsorption were determined  Clay mineral stabilized chitosan in acidic medium  Composites removed anionic dyes in acidic medium  The dye adsorption mechanism was understood via spectroscopic characterization Abstract Chitosan/KSF-montmorillonite beads were prepared at 1-25% w/w ratios by reaction of the biopolymer and the clay mineral with sodium tripolyphosphate as crosslinking agent. The synthesized biocomposites were characterized by CHN elemental analysis, XRD (X-Ray diffraction), FTIR (Fourier transform infrared), 13 C NMR (Nuclear magnetic resonance), thermogravimetry, SEM (Scanning electron microscopy), TEM (Transmission electron microscopy) and measurements of point of zero charge (pH pzc ).

The results suggested that hybrid materials were properly obtained and also their properties were improved compared to the pristine chitosan and montmorillonite.

Further, the beads were evaluated in the adsorption of Remazol Blue under batch operations obtained for different conditions at pH (2-8), contact time (0-660 min) and dye concentration (100-1600 mg L -1 ). Kinetic and equilibrium sorption parameters were evaluated using pseudo-first and pseudo-second rate models and Langmuir and Freundlich equations. The results were promising and suggested that the sorbents can be applied for removing anionic dyes from wastewater even in acidic medium.

Introduction

The presence of dyes in water is a serious environmental problem due their carcinogenic and mutagenic effects on aquatic life and also on human life [START_REF] Wang | Toxicity evaluation of reactive dyestuffs, auxiliaries and selected effluents in textile finishing industry to luminescent bacteria Vibrio fischeri[END_REF][START_REF] Sponza | Toxicity studies in a chemical dye production industry in Turkey[END_REF].

Dyes generally have a synthetic origin and a complex structure which ensures physico-chemical, thermal and optical stability due the presence of aromatic groups [START_REF] Langhals | Color Chemistry. Synthesis, Properties and Applications of Organic Dyes and Pigments. 3rd revised edition[END_REF].

Therefore, the removal of dyes from waste effluents is widely focused and the various techniques studied, among them adsorption has been shown to be one of the most applicable technologies [START_REF] Khajeh | Nanoadsorbents: Classification, Preparation, and Applications (with Emphasis on Aqueous Media)[END_REF][START_REF] Rouquerol | Adsorption by clays, pillared clays, zeolites and aluminophosphates[END_REF]. Various adsorbents such as clay minerals [START_REF] Zhu | Organo-clays as sorbents of hydrophobic organic contaminants: Sorptive characteristics and approaches to enhancing sorption capacity[END_REF][START_REF] Yuan | Clays and clay minerals for pollution control[END_REF][START_REF] Zhu | Adsorbents based on montmorillonite for contaminant removal from water: A review[END_REF] and biopolymers such as chitosan [START_REF] Ali | New Generation Adsorbents for Water Treatment[END_REF] have been used to remove different dyes from wastewater.

Chitosan (CS) is a natural polyaminosaccharide synthesized from the deacetylation of chitin, which is a polysaccharide consisting predominantly of unbranched chains of β-(1-4)2-acetoamido-2-deoxy-d-glucose [START_REF] Monteiro | The influence of chitosans with defined degrees of acetylation on the thermodynamic data for copper coordination[END_REF]. The biopolymer is widely used because of its special characteristics such as abundance in nature, nontoxicity, biocompatibility, biodegradability, hydrophilicity, antimicrobial activity and adsorption properties [START_REF] Wan | Adsorption of dyes and heavy metal ions by chitosan composites: A review[END_REF][START_REF] Hsu | Biocompatibility and antimicrobial evaluation of montmorillonite/chitosan nanocomposites[END_REF]. Active sites as amino and hydroxyl groups on CS make it a good candidate to interact with anionic dyes [START_REF] Kumar | Chitosan chemistry and pharmaceutical perspectives[END_REF]. Electrostatic interactions are reported to be involved in solid/liquid interface between protonated amino functions of the chitosan and the anionic groups of the dyes in solution. As consequence, this biopolymer has been widely studied in the removal of heavy metal cations and organopollutants from wastewater [START_REF] Wu | A review and experimental verification of using chitosan and its derivatives as adsorbents for selected heavy metals[END_REF][START_REF] Bekçi | Sorption of malachite green on chitosan bead[END_REF]. However, the use of CS in adsorption is pHdependent because it is soluble or form gel under pH below than 4 [START_REF] Pereira | Chitosan-montmorillonite biocomposite as an adsorbent for copper (II) cations from aqueous solutions[END_REF]. The chemical modification of chitosan and its hybrid derivatives, as crosslinking by tripolyphosphate (TPP) for instance, improve its chemical and mechanical properties [START_REF] Martins | Chitosan/TPP microparticles obtained by microemulsion method applied in controlled release of heparin[END_REF][START_REF] Ngah | Adsorption characterization of Pb(II) and Cu(II) ions onto chitosan-tripolyphosphate beads: Kinetic, equilibrium and thermodynamic studies[END_REF]. Therefore, the synthesis of based-chitosan composites has presented a new approach to modify various properties of this biopolymer, such as swelling/water adsorption [START_REF] Pongjanyakul | Investigation of novel alginate-magnesium aluminum silicate microcomposite films for modified-release tablets[END_REF] and mechanical/thermal performance [START_REF] Wang | Biopolymer chitosan/montmorillonite nanocomposites: Preparation and characterization[END_REF][START_REF] Wu | Biodegradable poly(lactic acid)/chitosan-modified montmorillonite nanocomposites: Preparation and characterization[END_REF]. In addition, resulting composites have a high rigidity, low specific gravity and a higher resistance to corrosion and oxidation [START_REF] Saravanan | Sorption studies on heavy metal removal using chitin/bentonite biocomposite[END_REF].

Various materials have been applied to prepare CS composites such as clay minerals because of their green and simple way of preparation and promising physical and chemical characteristics [START_REF] Darder | Chitosan-clay nanocomposites: application as electrochemical sensors[END_REF][START_REF] Günister | Synthesis and characterization of chitosan-MMT biocomposite systems[END_REF][START_REF] Wang | A novel biopolymer/rectorite nanocomposite with antimicrobial activity[END_REF][START_REF] Zhang | Preparation and characterization of chitosan nanocomposites with vermiculite of different modification[END_REF][START_REF] Nesic | Characterization of chitosan/montmorillonite membranes as adsorbents for Bezactiv Orange V-3R dye[END_REF]. Among clay minerals, montmorillonite (Mt) is the most used in the smectite group [START_REF] Jaber | Selectivities in Adsorption and Peptidic Condensation in the (Arginine and Glutamic Acid)/Montmorillonite Clay System[END_REF][START_REF] Jaber | Synthesis of Some Clay Minerals[END_REF]. It is constituted of a hydrous layered aluminum silicate with exchangeable cations and active hydroxyl groups. The nanometer thick layers of Mt consist of aluminum octahedron sheet sandwiched between two silicon tetrahedron sheets. Isomorphic substitution of Al 3+ with Mg 2+ in the octahedral sheets gives an overall negative charge, which is counterbalanced by cations (Na + , Ca 2+ , Mg 2+ and Fe n+ ) located in the interlayer space [START_REF] Jaber | Synthesis of Some Clay Minerals[END_REF][START_REF] Zhang | Synthesis of clay minerals[END_REF][START_REF] Adraa | Montmorillonite-Cysteine Composites for Heavy Metal Cation Complexation : a Combined Experimental and Theoretical Study[END_REF][START_REF] Jaber | Mercaptopropyl Al-Mg phyllosilicate: Synthesis and characterization by XRD, IR, and NMR[END_REF]. Clay minerals are abundant in nature and consequently low cost. Their properties such as the high cation exchange capacity, the possibility of interlayer expansion and the large specific area, may provide greater versatility in the interaction and intercalation of bulky molecules such as CS into Mt [START_REF] Monvisade | Chitosan intercalated montmorillonite: Preparation, characterization and cationic dye adsorption[END_REF][START_REF] Kittinaovarat | Chitosan/modified montmorillonite beads and adsorption Reactive Red 120[END_REF][START_REF] Jaber | Synthesis, characterization and applications of 2 : 1 phyllosilicates and organophyllosilicates: Contribution of fluoride to study the octahedral sheet[END_REF][START_REF] Jaber | New range of Al-Mg organoclays with tailored hydrophobicity: incorporation of fluoride as a local probe to study the octahedral character[END_REF][START_REF] Jaber | A new Al,Mg-organoclay[END_REF][START_REF] Gallego | Synthesis of new lamellar inorganic-organic talc-like hybrids[END_REF][START_REF] Jaber | Formation of organoclays by a one step synthesis[END_REF]. Although montmorillonite shows good interaction with cationic dyes [START_REF] Kittinaovarat | Chitosan/modified montmorillonite beads and adsorption Reactive Red 120[END_REF][START_REF] Wang | Removal of methylene blue from aqueous solution using chitosan-g-poly(acrylic acid)/montmorillonite superadsorbent nanocomposite[END_REF], it usually displays low affinity for anionic dyes because its negatively charged surface [START_REF] Fournier | Physicochemical characterization of lake pigments based on montmorillonite and carminic acid[END_REF][START_REF] Thangaraj | Fluorescent Quenching of Sulforhodamine dye over Graphene Oxide and Boron Nitride like Graphene[END_REF]. Therefore, the montmorillonite modification with polymers changes its physical and chemical properties and improves its affinity for anionic species [START_REF] Wan | Adsorption of dyes and heavy metal ions by chitosan composites: A review[END_REF][START_REF] Nesic | Characterization of chitosan/montmorillonite membranes as adsorbents for Bezactiv Orange V-3R dye[END_REF][START_REF] Monvisade | Chitosan intercalated montmorillonite: Preparation, characterization and cationic dye adsorption[END_REF]. However, a toxicological evaluation is always required to avoid potential human and environmental risks. Even though unmodified Mt has been reported to cause low cytotoxicity in many cell lines [START_REF] Maisanaba | Toxicological evaluation of clay minerals and derived nanocomposites: A review[END_REF], chitosan/TPP and chitosan/Mt composites has also been described to decrease the toxicity of Mt and to increase some CS properties, such as antimicrobial activity, better cell proliferation and biocompatibility [START_REF] Hsu | Biocompatibility and antimicrobial evaluation of montmorillonite/chitosan nanocomposites[END_REF]. Therefore, the CS-Mt composites are supposed to have very reduced toxicity.

In this perspective, the present investigation aimed to synthesize chitosanmontmorillonite biocomposites beads using TPP as crosslinked agent with Mt w/w ratios of 1, 5, 15 and 25% in order to obtain new adsorbents applicable for Remazol blue removal from aqueous solution. The Remazol blue RN (1-amino-9,10-dioxo-4-[3-(2-sulfonatooxietilsulfonil)aniline] anthracene-2-sulfonate, disodium) is an anionic dye widely applied in the textile industry, used as starting material to produce polymeric dyes [START_REF] Sathishkumar | Utilization of agro-industrial waste Jatropha curcas pods as an activated carbon for the adsorption of reactive dye Remazol Brilliant Blue R (RBBR)[END_REF] and is a toxic and recalcitrant organopollutants [START_REF] Mechichi | Remazol Brilliant Blue R decolourization by the laccase from Trametes trogii[END_REF][START_REF] Salleh | Cationic and anionic dye adsorption by agricultural solid wastes: A comprehensive review[END_REF]. The batch operations were evaluated at wide range of experimental conditions to investigate the effects of pH, contact time and the initial dye concentration on the sorption.

Materials and methods

Montmorillonite

Preparation of sodium KSF (KSF-Na)

Pristine montmorillonite-KSF was treated with a 0.1 mol dm -3 NaCl for 24 h at 323 K as described previously [START_REF] Gallego | Synthesis of new lamellar inorganic-organic talc-like hybrids[END_REF]. The cation exchange capacity of the KSF-Na was 60 cmol(+) Kg -1 [START_REF] Gallego | Synthesis of new lamellar inorganic-organic talc-like hybrids[END_REF].

Preparation of KSF-montmorillonite-chitosan beads

KSF-Na was suspended in 150 cm 3 of deionized water and the suspension was mechanically stirred for 24 h without heating. The w/w percentages of Na-KSF were 1, 5, 15 or 25% regarding the chitosan amount in the composite beads. To the clay suspension, 3.0 cm 3 of acetic acid and 3.0 g of CS were added and stirred for 24 h at room temperature. The gel was dropped slowly in the m/v 10% TPP solution. The KSF-CTS beads were matured for 24 h in the crosslinker solution, with distilled water up to neutral pH and dried at room temperature.

Dye removal

The influence of the pH on dye removal was investigated at the 2-8 pH range.

Samples of 50 mg of the beads were then suspended in 20.0 cm 3 dye solution at 500 mg L -1 and were maintained under agitation at 298 ± 1 K. The system was centrifuged and the residual dye concentration (C e ) was determined at 592 nm in a Shimadzu UV-Vis spectrometer, CBT-240 model. The amount of adsorbed dye was calculated by

Equation 1, (1)
where C i and C e are the initial and equilibrium dye concentrations (mg L -1 ), V is the solution volume (L) and m is the mass of the solid (mg).

For the sorption kinetic study, the same procedure was followed for times until 660 minutes at the optimal pH. For the equilibrium isotherm, initial dye concentrations varied from 100 to 1600 mg L -1 , under optimum pH and established time.

To ensure reproducible measurement results, the adsorption experiments were carried out in triplicate.

Characterizations

CHN elemental analyzes were performed at a Perkin-Elmer analyzer, PE Model 2400. XRD patterns of the powdered samples were performed using a Shimadzu XD3A diffractometer. FTIR spectra were collected with a Bomem MB-series spectrophotometer using KBr pellets. Thermogravimetric analysis was performed by a Shimadzu TGA-50H thermobalance. Surface morphology of samples was determined SEM in a JEOL JSTM-300 microscope. TEM was performed on a Philips CM 200 microscope operating at 200 kV. 13 C CP-MAS NMR spectra were obtained on a Bruker Avance 500 spectrometer operating at 125 MHz. The point of zero charge (pH PZC ) was determined after grinding the samples as previously described [START_REF] Silva | Adsorption of an industrial anionic dye by modified-KSF-montmorillonite: Evaluation of the kinetic, thermodynamic and equilibrium data[END_REF].

Results and Discussion

Characterization of chitosan-KSF hybrids

Composition, structural, thermal and morphologic characterizations

CHN results (Table 1) were used to estimate the content of clay mineral in the composites (Details see Table SM 1 15 and 25%) suggesting the clay mineral was totally incorporated in the beads. The quantities of organic groups (q N ) in the composites were estimated based on nitrogen content and the high obtained values for all composites were in agreement with the amount of chitosan in the beads. Table 1. CHN content for the pristine chitosan and their composites with KSF-Na, the quantity of organic groups on the solids (q N ) and estimated degrees of deacetylation (DD) for chitosan.

Sample

C H (%) N (%) q N (mmol g The IR spectrum of KSF-Na (Figure 1i-a), showed hydroxyl stretching bands located at 3634 cm -1 assigned to M-OH stretching (M = Mg, Al) of the octahedral sheets and at 3446 cm -1 and 1645 cm -1 attributed to the stretching and bending vibrations of the hydroxyl groups of interlayer water, respectively [START_REF] Madejová | FTIR spectroscopic characterization of thermally treated Cu2+, Cd2+, and Li+ montmorillonites[END_REF][START_REF] Tyagi | Determination of structural modification in acid activated montmorillonite clay by FT-IR spectroscopy[END_REF]. The bands between 1120 and 470 cm -1 were associated to Si-O-Si stretching and Si-O-Al deformation, respectively [START_REF] Tyagi | Determination of structural modification in acid activated montmorillonite clay by FT-IR spectroscopy[END_REF][START_REF] Li | An FTIR investigation of hexadecyltrimethylammonium intercalation into rectorite[END_REF]. The characteristic Si-O asymmetric and symmetric stretchings appeared at 1100 and 1040 cm -1 , respectively [START_REF] Zhang | Influence of four different types of organophilic clay on the morphology and thermal properties of polystyrene/clay nanocomposites prepared by using the γ-ray irradiation technique[END_REF][START_REF] Madejová | Comparative FT-IR study of structural modifications during acid treatment of dioctahedral smectites and hectorite[END_REF][START_REF] Paluszkiewicz | FTIR analysis of bentonite in moulding sands[END_REF], and the Si-O-Si deformation was detected at 525 cm -1 [START_REF] Tyagi | Determination of structural modification in acid activated montmorillonite clay by FT-IR spectroscopy[END_REF]. The bands between 916 (Al 2 OH) and 840 cm -1 (AlMgOH) corresponded to the octahedral layer of the clay mineral and reflected the partial substitution of Al octahedral by Mg [START_REF] Tyagi | Determination of structural modification in acid activated montmorillonite clay by FT-IR spectroscopy[END_REF][START_REF] Madejová | Comparative FT-IR study of structural modifications during acid treatment of dioctahedral smectites and hectorite[END_REF]. The band at 631 cm -1 was associated to quartz [START_REF] Tyagi | Determination of structural modification in acid activated montmorillonite clay by FT-IR spectroscopy[END_REF][START_REF] Paluszkiewicz | FTIR analysis of bentonite in moulding sands[END_REF].

The infrared spectrum of chitosan (Figure 1i-b) showed characteristic bands at 2920 and 2850 cm -1 assigned to C-H asymmetric and symmetric stretchings and a strong broad band around 3400 cm -1 assigned to the overlapping of O-H and N-H stretchings [START_REF] Domard | Preparation and characterization of fully deacetylated chitosan[END_REF][START_REF] Beil | Determination of the degree of N-acetylation (DA) of chitin and chitosan in the presence of water by first derivative ATR FTIR spectroscopy[END_REF][START_REF] Kasaai | A review of several reported procedures to determine the degree of N-acetylation for chitin and chitosan using infrared spectroscopy[END_REF]. The band around 1655 cm -1 was associated with the NH stretching and the C=O axial stretching from residual amide [START_REF] Domard | Preparation and characterization of fully deacetylated chitosan[END_REF][START_REF] Kasaai | A review of several reported procedures to determine the degree of N-acetylation for chitin and chitosan using infrared spectroscopy[END_REF]. Other characteristic bands were observed at: 1558 cm -1 , associated with NH 3 + deformation [START_REF] Wang | Biopolymer chitosan/montmorillonite nanocomposites: Preparation and characterization[END_REF][START_REF] Nesic | Characterization of chitosan/montmorillonite membranes as adsorbents for Bezactiv Orange V-3R dye[END_REF][START_REF] Paluszkiewicz | FT-IR study of montmorillonite-chitosan nanocomposite materials[END_REF], 1420 and 1320 cm -1 which were assigned to the C-N axial deformations from -CH 2 and amide III bands, respectively [START_REF] Beil | Determination of the degree of N-acetylation (DA) of chitin and chitosan in the presence of water by first derivative ATR FTIR spectroscopy[END_REF][START_REF] Kasaai | A review of several reported procedures to determine the degree of N-acetylation for chitin and chitosan using infrared spectroscopy[END_REF]; and at 1380 cm -1 which was attributed to CH-deformation of methyl group in N-acetyl [START_REF] Beil | Determination of the degree of N-acetylation (DA) of chitin and chitosan in the presence of water by first derivative ATR FTIR spectroscopy[END_REF][START_REF] Kasaai | A review of several reported procedures to determine the degree of N-acetylation for chitin and chitosan using infrared spectroscopy[END_REF]. Other bands were located at 1072 and 1029 cm -1 , due to glucopyranosidic ring stretchings [START_REF] Kasaai | A review of several reported procedures to determine the degree of N-acetylation for chitin and chitosan using infrared spectroscopy[END_REF][START_REF] López | A kinetic study on the thermal behaviour of chitosan[END_REF]; 1153 cm -1 , associated with betaglycosidic bond between carbons 1 and 4 and 890 cm -1 also associated with polysaccharide structure [START_REF] Monteiro | The influence of chitosans with defined degrees of acetylation on the thermodynamic data for copper coordination[END_REF][START_REF] Wan | Adsorption of dyes and heavy metal ions by chitosan composites: A review[END_REF][START_REF] Kittinaovarat | Chitosan/modified montmorillonite beads and adsorption Reactive Red 120[END_REF][START_REF] Laus | Adsorption and desorption of Cu(II), Cd(II) and Pb(II) ions using chitosan crosslinked with epichlorohydrintriphosphate as the adsorbent[END_REF]].

The characteristic bands of both KSF-Na and CS were observed in the spectra of beads (Figure 1i-c-f). However, the displacement of N-H deformation band at 1558 cm -1 for lower frequencies in the composites (1541 cm -1 ) suggested the interaction between the protonated amino groups of chitosan, TPP and the negative charge on the silicate layers [START_REF] Darder | Chitosan-clay nanocomposites: application as electrochemical sensors[END_REF][START_REF] Monvisade | Chitosan intercalated montmorillonite: Preparation, characterization and cationic dye adsorption[END_REF][START_REF] Tan | A novel method to prepare chitosan/montmorillonite nanocomposites in the presence of hydroxy-aluminum oligomeric cations[END_REF][START_REF] Celis | Montmorillonite-chitosan bionanocomposites as adsorbents of the herbicide clopyralid in aqueous solution and soil/water suspensions[END_REF].

Although there are several procedures to determine DD average for chitosan [START_REF] Kasaai | Determination of the degree of N-acetylation for chitin and chitosan by various NMR spectroscopy techniques: A review[END_REF][START_REF] Kasaai | Various Methods for Determination of the Degree of N-Acetylation of Chitin and Chitosan: A Review[END_REF][START_REF] Khan | Reporting degree of deacetylation values of chitosan: the influence of analytical methods[END_REF], IR spectroscopy is widely used [START_REF] Kasaai | A review of several reported procedures to determine the degree of N-acetylation for chitin and chitosan using infrared spectroscopy[END_REF][START_REF] Dimzon | Degree of deacetylation of chitosan by infrared spectroscopy and partial least squares[END_REF]. Therefore, Equation 6 [START_REF] Lima | A thermodynamic investigation on chitosan-divalent cation interactions[END_REF] was applied to determine DD for pure chitosan (CS):

DD = 97.67 -[26.486(A 1655 /A 3400 )] (6)
Thus, the DD value was 78 ± 1 for CS, in agreement with the estimated value by CHN.

KSF-Na diffraction patterns (Figure 1ii-a) presented characteristic reflection at 2θ = 6.29º, which was associated to the basal spacing (d 001 ) of 1.40 nm. Although the basal space has been reached a higher value than observed for montmorillonite around 1.2 nm [START_REF] Garrison | Molecular Modeling of Clay Mineral Structure and Surface Chemistry[END_REF][START_REF] Ait-Akbour | Adsorption of PolyCarboxylate Poly(ethylene glycol) (PCP) esters on Montmorillonite (Mmt): Effect of exchangeable cations (Na+, Mg2+ and Ca2+) and PCP molecular structure[END_REF][START_REF] Zhou | Modification of montmorillonite surfaces using a novel class of cationic gemini surfactants[END_REF], the same value was previously found [START_REF] Silva | Adsorption of an industrial anionic dye by modified-KSF-montmorillonite: Evaluation of the kinetic, thermodynamic and equilibrium data[END_REF][START_REF] Santos | Chitosan-Montmorillonite microspheres: A sustainable fertilizer delivery system[END_REF][START_REF] Liu | Adsorption of levofloxacin onto an iron-pillared montmorillonite (clay mineral): Kinetics, equilibrium and mechanism[END_REF], the higher value can be related to the acid treatment of raw KSF and the sodium exchange. Other characteristic peaks of montmorillonite occurred at 2θ ≅ 20 (d 020 = 0.45 nm) and 35 (d 006 = 0.26 nm).

Impurities of muscovite and quartz were also detected. formed by hydrogen bonds between the amino and hydroxyl groups of the CS chains [START_REF] Martins | Chitosan/TPP microparticles obtained by microemulsion method applied in controlled release of heparin[END_REF][START_REF] Samuels | Solid state characterization of the structure of chitosan films[END_REF]. The absence of reflection at 2θ ≅ 30 suggested CaCO 3 removal in the demineralization of the precursor chitin [START_REF] Younes | Chitin and Chitosan Preparation from Marine Sources[END_REF].

Q,M Q,M Q,M (f) (e) (d) (c) (b) (a)
The XRD patterns of the composites showed characteristic peaks of both KSF-Na and CS (Figure 1ii-c-f). For all prepared beads, the patterns presented broad peaks with low degree of crystallinity as observed in other studies [START_REF] Kittinaovarat | Chitosan/modified montmorillonite beads and adsorption Reactive Red 120[END_REF][START_REF] Wang | Removal of methylene blue from aqueous solution using chitosan-g-poly(acrylic acid)/montmorillonite superadsorbent nanocomposite[END_REF]. The reflection at 2θ = 6.29º (d 001 ) of montmorillonite almost disappeared, indicating highly disordered intercalated or probably exfoliated montmorillonite [START_REF] Wang | Biopolymer chitosan/montmorillonite nanocomposites: Preparation and characterization[END_REF][START_REF] Ray | 4 -Processing and Characterization[END_REF]. Also only broad reflections of CS were observed with maintenance of the peak at 2θ ≅ 20° in all samples.

The reaction between KSF-CS gels in the TPP solution involves the electrostatic interactions between protonated amino groups of the CS and the negative TPP, which resulted in the gelatinous beads. Thus, SEM studies only provide the information concerning the general morphology of the chitosan-clay composite. The surface chitosan-clay beads remained non-porous and became more irregular or roughened with the increasing of the KSF content in the composites (Figure 2).

Point of zero charge (pH PCZ )

KSF-Na surface is negatively charged in aqueous solution due to isomorphous substitutions of Al 3+ for Si 4+ in the tetrahedral layer and Mg 2+ for Al 3+ in the octahedral layer on montmorillonite [START_REF] Wan | Adsorption of dyes and heavy metal ions by chitosan composites: A review[END_REF][START_REF] Swartzen-Allen | Surface and colloid chemistry of clays[END_REF], (See Figure SM 1-a).

The isoelectric point of the CS was close to neutral (pH 6.5) (Figure SM 1-b). In acidic solution, the chitosan amino groups are protonated by H + by the solution, which increases the pH [START_REF] Rinaudc | Solubilization of Chitosan in Strong Acid Medium[END_REF]. The alkaline solution, the amino groups of chitosan are deprotonated, resulting in an excess of negative surface charge [START_REF] Kosmulski | Surface charging and points of zero charge[END_REF][START_REF] Neto | Preparation and evaluation of chitosan beads immobilized with Iron(III) for the removal of As(III) and As(V) from water[END_REF].

In the KSF-Na/CTS composites, the pH PCZ was more negative at lower pH as the KSF was added to beads, which indicated better chemical stability at acid medium. In fact, CS is soluble at pH lower than 4 forming gel, partially soluble at pH 4-6 and stable at pH higher than 6. CTS and KSF-CTS-1% were partially soluble at pH lower than 4 and stable at pH equal or higher than 4. The composites with 5-25% of KSF-Na became partially soluble at pH < 3 and stable at pH ≥ 3. These results indicated that the presence of the KSF-Na in the beads promoted a higher resistance to the dissolution of chitosan in acidic conditions.

Anionic dye removal study

Effect of pH solution

It is well-established that pH influences the adsorption of a given specie due to the variation of charges between the adsorbed species and the surface of the adsorbent [START_REF] Ngah | Adsorption characterization of Pb(II) and Cu(II) ions onto chitosan-tripolyphosphate beads: Kinetic, equilibrium and thermodynamic studies[END_REF][START_REF] Laus | Adsorption and desorption of Cu(II), Cd(II) and Pb(II) ions using chitosan crosslinked with epichlorohydrintriphosphate as the adsorbent[END_REF][START_REF] Kosmulski | Surface charging and points of zero charge[END_REF][START_REF] Baskaralingam | Modified hectorites and adsorption studies of a reactive dye[END_REF]. As UV-Vis spectrum is influenced by the solvent [START_REF] Montagner | Library of UV-Vis-NIR reflectance spectra of modern organic dyes from historic pattern-card coloured papers[END_REF][START_REF] Roca Jalil | Effect of the Al/clay ratio on the thiabendazol removal by aluminum pillared clays[END_REF], the behavior of dye at pH 3 to 8 was monitored (See Figure SM 2). Characteristics bands were observed in the UV-Vis spectra with partial overlap of two well-resolved bands between 450-700 nm, assigned to the π→π* and n→π* transitions of the dye [START_REF] Montagner | Library of UV-Vis-NIR reflectance spectra of modern organic dyes from historic pattern-card coloured papers[END_REF].

In order to check the optimum pH for adsorption, the process was monitored at the pH range of 2-8 (Figure 5). The RB sorption onto the beads was pH-dependent, which the maximum uptake of RB was achieved at pH 3. At lower pH, the primary amino groups of the CS and KSF-CTS composites beads become protonated (-NH 3 + )

and the surface of the beads became positively charged [START_REF] Roca Jalil | Effect of the Al/clay ratio on the thiabendazol removal by aluminum pillared clays[END_REF][START_REF] Ahmad | Equilibrium, Kinetics, and Thermodynamics of Remazol Brilliant Blue R Dye Adsorption onto Activated Carbon Prepared from Pinang Frond[END_REF][START_REF] Jiang | Adsorption of C.I. Reactive Blue 19 from aqueous solutions by porous particles of the grafted chitosan[END_REF]. Simultaneously, sulfonate groups of the dye (Dye-SO 3 Na) were dissociated in an aqueous solution (Dye-SO 3 -) [START_REF] Jiang | Adsorption of C.I. Reactive Blue 19 from aqueous solutions by porous particles of the grafted chitosan[END_REF]. Therefore, electrostatic interactions were favorable between the positively charged surface and the anionic RB dye at pH 3. This implies that strong electrostatic attraction is the dominant role for uptake of RB dye onto the given adsorbents, which can be favor the possible reuse of the adsorbents after RB desorption in alkaline medium.

Although the CS presented sorption capacity equivalent to the prepared composites, they proved to be inadequate to the sorption process when the pH is lower than 3.0 due to their partial dissolution, which results in loss of the adsorbent during the process. Moreover, no significant interactions occurred between KSF-Na and RB in the pH studied.

Effects of time on adsorption

The time isotherms are required for the preliminary study of the adsorptive process because they provide information about the kinetic for complete interaction to the surface of the adsorbent. Chitosan and composites were used with a 500 mg L -1 dye solution, under pre-determined pH previously obtained (Figure 6). Although time isotherms for RB dye reached equilibrium at 480 min for the prepared composites beads, it can be seen that at least 77% (160 mg g -1 ) of dye have been removed for all compositions at 240 min. Similar results were found for other based-chitosan materials [START_REF] Wan | Adsorption of dyes and heavy metal ions by chitosan composites: A review[END_REF][START_REF] Wang | Removal of methylene blue from aqueous solution using chitosan-g-poly(acrylic acid)/montmorillonite superadsorbent nanocomposite[END_REF][START_REF] Wang | Fast removal of methylene blue from aqueous solution by adsorption onto chitosan-g-poly (acrylic acid)/attapulgite composite[END_REF][START_REF] Auta | Chitosan-clay composite as highly effective and low-cost adsorbent for batch and fixed-bed adsorption of methylene blue[END_REF][START_REF] Mirmohseni | Chitosan hollow fibers as effective biosorbent toward dye: Preparation and modeling[END_REF], i.e. Mirmohseni, Seyed Dorraji, Figoli and Tasselli [START_REF] Mirmohseni | Chitosan hollow fibers as effective biosorbent toward dye: Preparation and modeling[END_REF] reached a RB dye removal efficiency around 80% at 240 min and pH 3.5 into chitosan hollow fibers biosorbents. The experimental data were fitted to pseudo-first order [START_REF] Lagergren | On the theory of so-called adsorption of materials[END_REF] and pseudo-second order models [START_REF] Ho | Pseudo-second order model for sorption processes[END_REF], as described the Equations 7 and 8,

  t k q q q e t e 303 . 2 log log 1    (7)
Where q e and q t are the amount of adsorbed dyes in equilibrium and at a time t and k 1 is the constant pseudo-first-order rate.

e e t q t q k q t   2 2 1 (8)
Where k 2 is the rate constant of pseudo-first order rate.

The results indicated that adsorption fitted well with the pseudo-second kinetic model with a high value of coefficient of determination (R 2 ) and a close agreement between experimental and estimated qe values (Figures SM 3, SM 4 and Table 2). It is known that the rate-limiting step for this model is the sorption mechanism involving chemisorption. Although passive diffusion agree with the time found to reach equilibrium in the isotherms, some other reaction mechanism may also been involved in the removal of the present dye. In addition, the rate constant (k 2 ) was affected by the KSF content in composites, which faster adsorption kinetics were observed for the composites with higher KSF quantity. Other chitosan-clay composites showed the same behavior for dyes removal [START_REF] Nesic | Characterization of chitosan/montmorillonite membranes as adsorbents for Bezactiv Orange V-3R dye[END_REF][START_REF] Peng | Chitosan-modified palygorskite: Preparation, characterization and reactive dye removal[END_REF][START_REF] Vanaamudan | Equilibrium, kinetics thermodynamic study on adsorption of reactive blue-21 and reactive red-141 by chitosan-organically modified nanoclay (Cloisite 30B) nano-bio composite[END_REF][START_REF] Liu | Adsorption of an anionic azo dye by crosslinked chitosan/bentonite composite[END_REF].

Table 2 -Kinetic parameters for the Remazol blue adsorption on chitosan samples and its composites with KSF-Na, according to the models of pseudo-first order and pseudosecond order at 500 mg L -1 , pH 3 and 298 ± 1 K.

Effect of the dye concentration

The equilibrium isotherms for RB are presented in Figure 7. Although KSF does not present any affinity to RB due to its negatively surface, the q e values were varied between 289 to 310 mg g -1 and suggest that the composition of the composites applied did not highly affect the adsorption capacity of the composites. According to literature dye adsorption on based chitosan/montmorillonite materials was influenced by the ratios of chitosan to clay mineral [START_REF] Darder | Chitosan-clay nanocomposites: application as electrochemical sensors[END_REF][START_REF] Zhang | Synthesis of clay minerals[END_REF][START_REF] Kittinaovarat | Chitosan/modified montmorillonite beads and adsorption Reactive Red 120[END_REF][START_REF] Jiang | Adsorption of C.I. Reactive Blue 19 from aqueous solutions by porous particles of the grafted chitosan[END_REF]. In the present study, based on CHN results, the lower quantity of chitosan in the beads was 3.83 mmol g -1 in the KSF-CTS-25%, which value is seven times higher than the maximum dye removal in the composite (310 mg g -1 or 0.5 mmol g -1 ).

Furthermore, the increase of the KSF amount in the prepared composites has contributed to obtain a low-cost adsorbent for anionic dye with a higher resistance in acidic solutions. Other point is the facility of separation of the beads after adsorption. The data were adjusted to Langmuir [START_REF] Langmuir | The adsorption of gases on plane surfaces of glass, mica and platinum[END_REF] and Freundlich models (Details see [START_REF] Freündlich | Über die Adsorption in Lösungen[END_REF]. The Langmuir (Equation 9) model was applied for describing the adsorption data:

Figure SM 5)
max max 1 bq q C q C e e e   (9)
where C e is the concentration of dye in equilibrium; q e is the amount of adsorbed dye at the solid/liquid interface and q max is the monolayer adsorption capacity of the adsorbent, both in terms of quantity of material by mass of function; b is the Langmuir adsorption constant associated with adsorption energy.

The Freundlich model [START_REF] Freündlich | Über die Adsorption in Lösungen[END_REF] is based on the heterogeneity or non-uniformity of sites on solid and the multilayer adsorption and is expressed as Equation 10, e F e C n k q ln 1 ln ln   [START_REF] Monteiro | The influence of chitosans with defined degrees of acetylation on the thermodynamic data for copper coordination[END_REF] where k F and 1/n are empirical parameters of the Freundlich model. k F is the adsorption constant, related to the binding energy, and 1/n is connected to the heterogeneity of the solid surface. The parameter n is a measure of the deviation from linearity of the adsorption and indicates the degree of non-linearity between the solution concentration and adsorption.

The resulted parameters are summarized in Table 3. Langmuir model showed higher R 2 than Freundlich model for the prepared composites and q max closed to q e values and exhibited an adsorption capacity around 310 mg g -1 for RB removal, suggesting the monolayer of the anionic dye on composites and that chemisorptions is the principal mechanism of interaction between the anionic dye and (-NH 3 + ) protonated surface of the beads.

Characterization of chitosan-KSF-anionic dye hybrids

The post-adsorption samples presented some modifications in FTIR spectra The UV-Vis spectrum for KSF-Na showed an absorption band centered on 260 nm [START_REF] Chen | Synthesis and characterization of TiO2-montmorillonites doped with vanadium and/or carbon and their application for the photodegradation of sulphorhodamine B under UV-vis irradiation[END_REF][START_REF] Alin | In situ characterization of organo-modified and unmodified montmorillonite aqueous suspensions by UV-visible spectroscopy[END_REF] and for chitosan two bands, a broad one at 300 nm and a sharp one at 400 nm (Figure SM 7) [START_REF] Iturriaga | Active naringinchitosan films: Impact of UV irradiation[END_REF][START_REF] Andrady | Spectral sensitivity of chitosan photodegradation[END_REF]. The spectra of the chitosan-KSF/dye hybrids presented new broad band at around 600 nm associated probably with the bathochromic effect of the original absorption of the free dye at 543 nm.

Conclusions

Beads chitosan/-montmorillonite-KSF at 1, 5, 15 and 25% w/w proportions were prepared using sodium tripolyphosphate as crosslinker. XRD, FTIR and SEM results highlighted the interactions between the biopolymer and the clay mineral leading to a hybrid composite with better adsorption properties for the anionic dye. TEM showed co-existence of separated phases and exfoliated structures. Measurements of the zero charge potential (pH ZCP ) suggested that isoelectric point shifted to lower pH with the increasing of montmorillonite content, which extends the application at more acidic medium.

The optimal conditions for the removal of Remazol Blue-RN from aqueous solution were pH 3 during 480 min. The processes were followed pseudo-second order kinetics, with the rate constant dependent on montmorillonite content in the composites.

The equilibrium isotherms indicated that the Remazol blue sorption capacity was substantially independent of clay mineral content in the beads.

Thus, the prepared composites showed better adsorption properties than pure montmorillonite with maximum dye removal of 310 mg g -1 for the composite prepared with 25% of montmorillonite.

This study demonstrated that chitosan/KSF beads can be applied as versatile bioadsorbent for anionic dyes removal from solution in a wide range of pH and especially in acidic medium where the chitosan is soluble.

Scheme 1 .

 1 Scheme 1. Chemical structure of Remazol Blue RN.

, ES 1 ,

 1 Equations 2 to 5) and DD average once the crosslinking with TPP does not affect the quantity of acetylated and deacetylated units of the chitosan. In other words, the DD should not be affected in the various composites. The estimated values varied around 78 ± 2 in the composites. The pristine KSF-Na showed C<0.3% and N<0.05%. The carbon and nitrogen contents decreased proportionally with the increasing of KSF-Na in the composite. The experimental quantities of the clay mineral were in agreement with the proposal compositions(1, 5,

Figure 1 .

 1 Figure 1. i) Infrared spectra and ii) for X-ray diffractograms for (a) KSF-Na, (b) CS,

Figure 2 .Figure 3 .

 23 Figure 2. SEM images for (a) KSF-CTS-1%, (b) KSF-CTS-5%, (c) KSF-CTS-15% and

Figure 4 .

 4 Figure 4. 13 C NMR spectra for (a) CS, (b) KSF-CTS-1%, (c) KSF-CTS-5%, (d) KSF-

Figure 5 .

 5 Figure 5. pH effect over Remazol blue adsorption for (a) CS, (b) KSF-CTS-1%, (c) KSF-CTS-5%, (d) CTS-KSF-15% and (e) KSF-CTS-25% samples at 500 mg L -1 , 24 h

Figure 5 .

 5 Figure 5. pH effect over Remazol blue adsorption for (a) CS, (b) KSF-CTS-1%, (c) KSF-CTS-5%, (d) CTS-KSF-15% and (e) KSF-CTS-25% samples at 500 mg L -1 , 24 h

Figure 6 .

 6 Figure 6. Time effect over Remazol blue removal for (a) CS, (b) KSF-CTS-1%, (c) KSF-CTS-5%, (d) KSF -CTS-15% and (e) CTS-KSF-25% at 500 mg L -1 , pH 3 and 298

Figure 7 .

 7 Figure 7. Effect of the dye initial concentration on Remazol blue sorption for (a) CS,

(Figure SM 6 -

 6 i). The strong and broader O-H and N-H groups overlapped stretching bands, originally at 3425 cm -1 in pristine composites, became sharper and shifted to 3437 cm -1 due to the presence of these groups from the RB sorption. The characteristic bands for the dye (1624, 1591, 1529 and 1408 cm -1 )[START_REF] Hough | Trends in Solar Energy Research[END_REF] can overlapped with the composites structural bands. The band corresponding to ⱱs (Dye-SO 3 ) appears at 1192 cm -1 . Moreover, the band at 1558 cm -1 (in CS and KSF-CTS) assigned to -NH 3 + shifted to 1578 cm -1 for KSF-CTS-RB composites indicating a possible interaction between the -NH 3 + groups and the basic sites of the dye.The XRD patterns of the samples after RB adsorption (Figure SM 6-ii) showed similar reflections than those observed for montmorillonite and chitosan, suggesting the dye accommodation on surface of the beads.

Table 3 -

 3 Adsorption isotherms parameters for the Remazol blue on chitosan samples and their composites with KSF-Na, according to the models ofLangmuir and

	Freundlich at 100-1600 mg L -1 , pH 3 and 298 ± 1 K.
				Sample
	Model
		CS	KSF-CTS-1% KSF-CTS-5% KSF-CTS-15% KSF-CTS-25%
	Langmuir
	q e (mg g -1 )	308.34 ± 0.487 289.69 ± 0.890 306.18 ± 0.683 309.18 ± 0.493 310.73 ± 0.556
	q max (mg g -1 )	300.30	290.70	334.78	308.64	308.64
	b (dm g -1 )	0.2502	0.0606	0.1544	0.1192	0.1528
	R 2	0.9994	0.9979	0.9981	0.9979	0.9995
	Freündlich
	K f	7.5254	6.5914	7.7302	6.7924	6.3934
	n	6.8451	4.9848	7.7302	4.3630	3.8319
	R 2	0.4827	0.9522	0.7237	0.4399	0.6230

Table 4 )

 4 and demonstrated that all prepared composites in the present work had a good efficiency to RB uptake from wastewater.

Table 4 .

 4 Dye removal capacities of the chitosan-montmorillonite composites reported in this paper with similar adsorbents.

	Adsorbent	Dye	Concentration of dye (mg L -1 )	Adsorption capacity (mg g -1 )
	Chitosan	Remazol Blue RN	1700	308*
	Chitosan-modified palygorskite	Reactive Red 120	140	71 94
	Chitosan-Organically modified Nanoclay (Cloisite30B)	Reactive dyes RR-141	1000	443 95
	Chitosan-Organically modified Nanoclay (Cloisite30B)	Reactive dyes RB-21	1000	490 95
	Chitosan/Bentonite cross-linked with glutaraldehyde	Amido Black 10B	550	324 96
	Chitosan/MMT membrane	Bezactiv Orange	80	279 23
	Chitosan/montmorillonite-1% cross-linked with TPP	Remazol Blue RN	1700	290*
	Chitosan/montmorillonite-5% cross-linked with TPP	Remazol Blue RN	1700	306*
	Chitosan/montmorillonite-15% cross-linked with TPP	Remazol Blue RN	1700	309*
	Chitosan/montmorillonite-25% cross-linked with TPP	Remazol Blue RN	1700	310*
	*Present work

Acknowledgments

KSS was supported by a pos-doc grant from CNPq. CNPq is acknowledged for research fellowships to MGF and IMGS. CAPES/COFECUB (835/2015) for financial support.

SUPPLEMENTARY MATERIAL

Table SM1. Quantities of chitosan and KSF-Na on the beads based on CHN data


F Bico

Melo

B Darbois Texier

A Ibarra

F Vivanco

J Bico

F Melo

Friction of a sphere rolling down a granular slope

Keywords: 45.70.-n -Granular systems, classical mechanics of, 81.05.Rm -granular materials 07.79.Sp -friction force

Friction of a sphere rolling down a granular slope

Introduction. -The scientific questions associated with rolling friction on granular beds are crucial for understanding motorized locomotion on sandy environments, the evolution of talus [START_REF] De Blasio | Dynamics of grains falling on a sloping granular medium: application to the evolution of a talus[END_REF], ballistics of spherical asteroid landers [START_REF] Van Wal | Rolling resistance of a spherical pod on a granular bed[END_REF] and even the rolling of heavy dung balls in soft soils by beetles [START_REF] Bartholomew | Endothermy in african dung beetles during flight, ball making, and ball rolling[END_REF]. This phenomenon has benefited from a revival of interest after the incident of the rover Curiosity, developed by the NASA for Martian exploration, that got stuck for two weeks in a sand dune in 2005. Moreover, the need of robust technological solutions aimed at exploring deformable grounds has motivated the development of instrumented spheres that are propelled by winds similarly as tumbleweeds [START_REF] Antol | Low cost mars surface exploration: The mars tumbleweed[END_REF][START_REF] Wilson | Design, analysis, and testing of mars tumbleweed rover concepts[END_REF]. This solution has been also adapted to polar explorations [START_REF] Behar | Nasa/jpl tumbleweed polar rover[END_REF] and reforestation of arid zones. However, the success of these instrumented spheres in the exploration of complex environments requires a deep knowledge of the interaction between the non-cohesive substrates and the sphere. Historically, the problem of rolling friction on a deformable bed has been addressed by Bekker from the point of view of soils mechanics [START_REF] Bekker | Theory of land locomotion[END_REF][START_REF] Bekker | Off-the-road locomotion[END_REF]. Bekker expressed the drag of a wheel on a flat granular substrate as a function of its geometry, load and empirical soil constants deduced from soil tests. Although this approach has proven effective in practical situations, physical understanding is still missing. More recently, Van Wal et al. studied experimentally and numerically the rolling resistance of a light sphere on a horizontal granular bed [START_REF] Van Wal | Rolling resistance of a spherical pod on a granular bed[END_REF]. They developed a theoretical model based on the micro-collisions between the sphere and the ground to account for the rolling friction that partially agrees with their experimental observations. Beyond flat substrates, the case of sandy slopes has been studied by Crassous et al. in the context of ants captured inside antlions sand traps [START_REF] Crassous | Pressuredependent friction on granular slopes close to avalanche[END_REF]. They show that the friction of a slider on a granular bed exhibits large deviations from the Amontons-Coulomb laws at small applied pressure. Predicting the frictional behaviors on inclined granular slopes thus require new physical laws that remain to be found. The problem of the rolling resistance of dense spheres on an inclined sand bed has been addressed by De Blasio and Saeter [START_REF] De Blasio | Rolling friction on a granular medium[END_REF]. It was observed that the penetrating depth and the rolling friction of the spheres both increase with their density. Since this study was limited to a single slope angle, more experimental efforts are required to achieve a complete physical description of the problem. Here, we investigate systematically the dynamics of spheres rolling down an inclined surface of a non cohesive granular sample. The first section describes the experimental setup and the inferred sphere dynamics. The drag of the sphere is estimated as a function of the slope angle and the main physical parameters. Simultaneously, the deformation of the granular bed due to the passage of the sphere is characterized through the deflection of a laser sheet. In the second part, a simple mechanical model for the rolling resistance on a granular bed is proposed and compared to our experiments. This model accounts for a gradual increase of the asymmetry of the stress distribution acting on the area of the sphere in contact with the bead as the slope is increased. Finally, the different behaviors of the sphere on a granular bed are discussed and summarized in a phase diagram constructed with the slope angle and the ratio of the densities of the sphere and of the grains.

Experiments. -The experiment consists in releasing a homogeneous sphere of radius R and density ρ s over a flat granular bed made of glass particles of diameter d = 0.5 mm (polydispersity of 15%) and tilted by an angle α from the horizontal direction (Fig. 1). Spheres are made of expanded polystyrene foam and have a smooth surface in comparison to the grain size. The granular bed is prepared by pouring grains on a horizontal substrate covered with sandpaper and increasing its slope up to an angle α a = 23 • at which an avalanche starts. Once the avalanche stops, a granular layer of constant thickness h 10 mm remains and covers the substrate. According to the study of Pouliquen & Forterre [START_REF] Pouliquen | Friction law for dense granular flows: application to the motion of a mass down a rough inclined plane[END_REF], the height of the granular layer, h depends solely on the slope angle, which allows for the preparation of layers of well defined thickness and good homogeneity. We verified that increasing the thickness of the granular layer does not affect the dynamics of the sphere, indicating that our experiments correspond to the limit of deep granular beds (h > 4 mm). The bulk density of this preparation has been measured to be ρ g = 1540 kg/m 3 , corresponding to a packing fraction of φ = 0.58. Thereafter, the slope of the system is reduced to the desired angle α at which the experiment is conducted. The sphere is released at null height above the granular bed without any initial velocity or initial spin. Within the range of sphere radius and density, selected for the experiments, the sphere does not sink entirely into the granular layer and rolls down the incline. The dynamics of the sphere along the slope is captured at rate of 30 fps through a CCD camera located one meter above the granular layer, preventing optical distortions. The sphere motion is extracted from movies analysis through a MAT-LAB code based on image correlation.

The typical dynamics of a sphere rolling on the granular slope [Fig. 2] reveals that the falling acceleration increases with α and that data are well described by a quadratic law x = a t 2 /2 [solid lines in Fig. 2(a)], where a is the acceleration of the sphere that is approximatly constant during a fall. Similar experiments were reproduced with spheres of different diameters and densities. In all cases, good fit of experimental data was achieved and the mean acceleration of the sphere, a, was deduced.

Figure 2(b) indicates that a minimal slope angle α c is required for the sphere to initiate spontaneously its motion. For α > α c , the acceleration of the sphere increases with α and shows vanishing dependency on R [blue symbols in Fig. 2(b)]. Increasing the sphere density ρ s results in increasing α c and reducing the sphere acceleration a. For a given angle, the measured acceleration is significantly lower than the one expected for the classical problem of a sphere rolling without sliding [a = 5/7 g sin α, cf. dark solid line and white dots in Fig. 2(b)]. Thereafter, we investigate the rotation of the sphere while rolling down the granular slope by means of black painted stripes on the sphere surface. Tracking these stripes provides the sphere rotation angle θ recovered while it travels over a given distance x. We observe that the ratio Rθ/x is equal to unity within the experimental accuracy of this study [Fig. 3(a)]. We thus conclude that the sphere does not significantly slide. The difference in sphere acceleration compared to the case of a flat and rigid substrate is then attributed to the friction of the sphere due to the deformation of the granular bed. Introducing the friction force F d experienced by the sphere along the x-direction, the equation of motion yields,

ẍ g = sin α - F d M g , (1)
where M = 4πρ s R 3 /3 and g is the gravitational acceleration. Considering that the sphere acceleration ẍ is constant and equals to the fitted value a, Eq. ( 1) leads to an estimate of the granular friction F d . The rolling of the sphere on the granular bed pushes away excess granular material which results in the formation of a characteristic groove as discussed by Crassous et al. [START_REF] Crassous | Pressuredependent friction on granular slopes close to avalanche[END_REF]. In order to characterize the features of the groove, we follow the profile of the granular bed with a laser sheet that intercepts the surface with an oblique incidence, as depicted in Fig. 1. A top view of the groove boundary is shown in Fig. 4. We first observe that for a given ex- x periment, the depth of the groove δ and its width w, both remain unchanged along the path of the sphere. These results can be attributed to the fact that for α < 22 • the sphere displacement does not trigger avalanches but only a local flow of grains. Conversely, for 22 • < α < 23 • , the deposition of the sphere induces global avalanches, a regime that lays out of the scope of this study. We then investigate the dependence of δ and w on the radius of the sphere for a given density [Fig. 5 Model.

(mm) α = 10 • α = 12 • α = 14 • α = 16 • α = 18 • α = 20 • α = 22 • (b)
-We first focus on the features of a groove produced by a sphere that is gently deposited on a flat and horizontal granular surface. Uehara et al. investigated the depth and the diameter of a crater formed by a sphere impacting a granular medium at low speed [START_REF] Uehara | Low-speed impact craters in loose granular media[END_REF]. In the limit of sphere released just above the medium, the and the width of the crater read,

δ = A µ 3/2 g ρ s ρ g 3/4 R (2)
w = B µ 7/8 g ρ s ρ g 7/16 R (3)
where µ g is the friction coefficient between grains (for the glass beads used in our experiments, µ g = 0.40, as deduced from the angle of repose of a granular pile) and A, B are numerical coefficients. Adjusting Eqs. ( 2) and ( 3) to the experimental data, we find that the best fits are [START_REF] Uehara | Low-speed impact craters in loose granular media[END_REF]. We attribute this difference to the fact that the contact area between the rolling sphere and the grains is lower than the area of a sphere at rest.

Hereafter, we develop an empirical model for the granular friction experienced by the sphere. While rolling downwards, the sphere is in contact with grains on a zone characterized by the angle θ m given by cos θ m = 1 -δ/R [Fig. 7]. In a continuum approach, normal and tangential stresses are distributed along the contact zone. These stress distributions can be reduced to two punctual forces along the radial and orthoradial direction, denoted f r and f θ respectively. The orientation of f r and f θ relatively to the slope defines the angle θ c [Fig. 7]. The static equilibrium of the sphere thus corresponds to θ c = α, f r = -M g and f θ = 0. For the case of rolling motion, we propose a heuristic expression of θ c as a function of θ m and α, θ c = θ m 1 -e -c α/θm . This expression includes a progressive loss of symmetry of the stress distribution for increasing α above α c and the saturation of θ c to θ m for large enough α. Indeed, the center of mass of the stress distribution cannot exceed the size of the contact zone. Note that the proposed expression only includes a single adjustable parameter denoted c. Following this expression of θ c , the critical slope angle α c for which the sphere starts rolling respects α c /θ m = 1 -e -c αc/θm . Along the y direction, f r and f θ projections compensate the weight of the sphere, leading to,

-f r cos θ c + f θ sin θ c = M g cos α (4)
Along the x direction, the projections of f r and f θ define the friction force F d experienced by the sphere,

F d = -f r sin θ c -f θ cos θ c (5)
Finally, the torque equilibrium toward the center of the sphere yields:

J θ = R f θ (6)
where J = 2/5M R 2 is the moment of inertia of a homogeneous sphere relatively to the sphere centered axis. Measurements of the sphere rotation indicate that the nonsliding approximation holds in our experiments [Fig. 3(a)]. Thus, the non-sliding relation, R θ = -ẋ, combined with Eqs. ( 1) and ( 6) leads to,

f θ = - J M R 2 (M g sin α -F d ) , (7)
which along with Eqs. (4), (5) yields to the reduced friction force

F d M g = cos α 2 5 tan α + sin θ c 2 5 + cos θ c (8)
where θ c = θ m 1 -e -c α/θm , cos θ m = 1 -δ/R, and δ is given by Eq. ( 2). Equation [START_REF] Bekker | Off-the-road locomotion[END_REF] predicts that F d /M g increases with the slope angle α and the density ratio ρ s /ρ g and does not depend on the sphere radius R. These trends are in qualitative agreement with the observations made previously [Fig. 3(b)]. Thereafter, the experimental data have been adjusted to the theoretical predictions of the normalized friction force F d /M g deduced from Eq. ( 8).

The best fit of the data is provided for c = 1.5 as represented by blue and yellow solid lines in Fig. 3(b). The good agreement with experiments justifies a posteriori the assumption made on the dependance of θ c with θ m and α. Indeed, the experimental data could not be captured by an empirical law where θ c depends solely on θ m . Moreover, the fact that a same value of the coefficient c allows to fit different density ratios ρ s /ρ g suggests that c does not depend on the properties of the granular bed. In addition, we performed similar experiments on a sand bed (ρ g = 1549 kg/m 3 and µ g = 0.65) and a good agreement between the theory and the experiments has been obtained when the physical properties of the sand are considered.

Finally, the fit of the predictions with experimental data provides a first characterization of the distribution of contact forces and its dependency with the properties of the granular medium and the slope angle.

Discussion. -Following Eq. ( 8), we deduce that the dynamics of the rolling sphere on a deformable granular slope is equivalent to the dynamics of a sphere rolling over a non deformable incline with a friction coefficient µ eff = 2 5 tan α + sin θ c / 2 5 + cos θ c . This observation justifies De Blasio and Saeter assumption of interpreting the deceleration of a sphere on a granular slope in term of effective friction coefficients [START_REF] De Blasio | Rolling friction on a granular medium[END_REF]. In their experiments, an increase of the effective friction coefficient with the slope angle was observed, a fact which is in agreement with the predictions of Eq. ( 8). However, a direct comparison with our model is not possible because the internal friction coefficient of their grains was not provided. Note that the expression of µ eff depends on the properties of the granular medium through its bulk density and its internal friction coefficient. The larger the values of ρ g and µ g , the lower the normalized granular friction F d /M g. We thus expect the compaction of the granular medium to reduce the friction of the rolling sphere. This fact was verified experimentally by tapping several times the granular preparation prior to the release of the sphere.

Thereafter, we apply our results to estimate the force required for beetles to push dung balls over sandy surfaces. Considering ρ g ∼ 1500 kg/m 3 , µ g ∼ 0.7 and assuming that the dung ball density is ρ b ∼ 800 kg/m 3 , the force required to roll the dung ball on a flat and sandy ground is equal to about a fifth of the ball weight. For a dung ball that is 50 times the weight of the insect as measured by Bartholomew and Heinrich [3], a beetle must be able to generate forces about 10 times its own weight in order to move the dung ball. The maximal pushing force of beetles belonging to the Scarabaeidae family, has been reported by Evans and Forsythe to be about 250 times their own weight [START_REF] Evans | A comparison of adaptations to running, pushing and burrowing in some adult coleoptera: especially carabidae[END_REF], thus ball rolling is not the hardest tasks that dung beetles can complete.

In this work, we only considered the case of high con-trast in density where the rolling sphere induces a surrounding granular flow. However, our predictions do not hold for low density ratio where the penetration depth δ becomes lower than the grain size (δ < d), a fact that occurs for ρ s /ρ g < 11 µ 2 g (d/R) 4/3 according to Eq. (2). In this limit, the granular bed acts as a non-deformable rough substrate for the rolling sphere. Such a situation has been investigated by Tardivel et al. with balls decelerating on a flat granular bed [START_REF] Van Wal | Rolling resistance of a spherical pod on a granular bed[END_REF]. In their experiments, they used a basket ball of radius R = 119 mm and density ρ s = 89 kg/m 3 and a medicine ball of radius R = 113 mm and density ρ s = 447 kg/m 3 . The granular bed is made of grains of size ranging from 6 to 12 mm with an internal friction coefficient µ g = 0.87. In such conditions, the penetration depth δ predicted by Eq.( 2) is less than 2 mm, thus lower than a third of the minimal grain size. In order to approach this limit, Tardivel et al. developed a model based on a succession of micro-collisions on a rigid surface with small asperities. Such a model yields an expression for the effective friction coefficient µ eff = M R 2 /J d/2R which predicts the experimental value within a precision of 25 %.

Finally, the different behaviors of a sphere placed on a granular slope are summarized in Fig. 8. Stable granular slopes are encountered for angles smaller than the avalanche angle α a (below the vertical dashed line). For small density ratio ρ s /ρ g < 11 µ 2 g (d/R) 4/3 , the granular bed remains undeformed and behaves like a rough surface as described by Tardivel et al. [START_REF] Van Wal | Rolling resistance of a spherical pod on a granular bed[END_REF] (below the horizontal dotted line). At small slope angles α < α c or large density ratio, the sphere is in static equilibrium in the crater formed by its release (above the solid line which represents the solution of the relation α c /θ m = 1 -e -c αc/θm that leads to α c 0.61 θ m for c = 1.5). Between these three limits, the sphere rolls down while deforming the granular bed as considered in the present study. Conclusion. -This study investigated the dynamics of a sphere rolling down a granular slope. While moving down, the sphere displaces grains and leaves a track, inducing a granular friction. Remarkably, the sphere follows an uniformly accelerated motion within the range of experimental parameters explored in this work. We estimated the steady friction force from the dynamics of the sphere and showed that it increases with the sphere density and the slope angle. The depth and the width of the groove created by the sphere were characterized experimentally and approached by the laws established in the case of a sphere impacting a granular bed at zero falling height. Thereafter, a mechanical model of the rolling sphere on a granular slope has been developed. Although this model is based on a heuristic description of the distribution of the contact forces between the rolling sphere and the granular bed, the comparison with experiments is fairly good. Finally, this model provides crucial informations on the distribution of contact forces and its dependency on the slope angle. * * *

Fig. 1 :

 1 Fig. 1: Sketch of the experimental setup.

Figure 3 (

 3 b) shows F d /M g as a function of α for different sphere diameters and densities. As a general trend, F d increases both with α and ρ s and is larger than in the case of pure rolling on a flat and rigid substrate where F d = 2/7 M g sin α [black solid line and white dots in Fig.3(b)].

α

  ( • ) a (m/s 2 ) a = 5/7 g s in α Rigid s ubs t r at e R = 11 mm, ρ s = 22 kg/m 3 R = 16 mm, ρ s = 22 kg/m 3 R = 24 mm, ρ s = 22 kg/m 3 R = 22 mm, ρ s = 57 kg/m 3

Fig. 2 :

 2 Fig. 2: (a) Position of a sphere as a function of time for different slope angles α. The sphere diameter is R = 16 mm and its density ρs = 22 kg/m 3 . Grains are made of glass beads of mean diameter d = 0.5 mm (the polydispersity is about 15%) and the bulk density of the bed is ρg = 1540 kg/m 3 . The granular layer is prepared by bringing the system to α = 23 • and waiting for the layer to reach the equilibrium thickness of h 10 mm. Solid lines correspond to quadratic fits of the form x = at 2 /2.(b) Acceleration of the sphere as a function of the slope angle α for different sphere diameters (R = 11, 16, 24, 22 mm for colored dots, triangles, squares and diamonds respectively) and densities (ρs = 22, 57 kg/m 3 for blue and yellow symbols). The dark solid line represents the classical theory for rolling without sliding for a homogeneous sphere (a = 5/7 g sin α). The white dots correspond to the case of a sphere rolling on a rigid incline.

R

  (a)-(b)]. One notices that the groove depth and width both increase proportionally to R. Furthermore, the detailed study of the geometry of the groove as a function of the slope α for a given sphere = 11 mm, ρ s = 22 kg/m 3 R = 16 mm, ρ s = 22 kg/m 3 R = 24 mm, ρ s = 22 kg/m 3 R = 22 mm, ρ s = 57 kg/m 3 M g s in α 2/7 M g s in α Rigid s ubs t r at e

Fig. 3 :Fig. 4 :Fig. 5 : 2 g

 3452 Fig. 3: (a) Rolling ratio Rθ/x as a function of the slope α for a sphere of diameter R = 24 mm and density ρs = 22 kg/m 3 over a granular bed made of glass beads. The black solid line indicates the classical case of pure rolling. (b) Granular reduced friction F d /M g deduced from experimental accelerations and Eq. (1). Symbols have the same signification as in Fig. 2(b). The blue and yellow solid lines show the prediction of Eq. (8) for c = 1.5. Gray regions show the model predictions for c = 1.5 ± 0.2. The dark solid line represents the tangential force in the classical case of a pure-rolling motion on a rigid substrate (F d = 2/7M g sin α).

Fig. 6 :

 6 Fig. 6: Depth (a) and width (b) of the groove left after the passage of a sphere as a function of the slope angle α for a sphere of radius R = 22 mm and density ρs = 57 kg/m 3 . Squares correspond to experimental data and black solid lines to the predictions of Eqs. (2) and (3) for ρs/ρg = 57/1540, µg = 0.40, A = 0.16 and B = 1.94.

Fig. 7 :

 7 Fig. 7: Notations used in the model.

Fig. 8 :

 8 Fig.8: Diagram of the different regimes encountered by a sphere placed on a granular bed made of glass beads as a function of its angle α and the density ratio between the sphere and the grains ρs/ρg. Crosses and dots indicate respectively the non-rolling and rolling situations observed experimentally for spheres of respective density ρs = 21, 57 and 171 kg/m 3 .

The authors acknowledge Fivos Pham for preliminary experiments. F.M. acknowledges support from Proyecto Fortalecimiento USACH USA1799-MF072215.


Mathis Fleury

Christian Barillot

Marsel Mano

Elise Bannier

Pierre Maurel

March

Automated Electrodes Detection during simultaneous EEG/fMRI

The coupling of Electroencephalography (EEG) and functional magnetic resonance imaging (fMRI) enables the measure of brain activity at high spatial and temporal resolution. The localisation of EEG sources depends on several parameters including the knowledge of the position of the electrodes on the scalp. An accurate knowledge about this information is important for source reconstruction. Currently, when acquiring EEG and fMRI together, the position of the electrodes is generally estimated according to fiducial points by using a template. In the context of simultaneous EEG/fMRI acquisition, a natural idea is to use magnetic resonance (MR) images to localise EEG electrodes. However, most MR compatible electrodes are built to be almost invisible on MR Images. Taking advantage of a recently proposed Ultra short Echo Time (UTE) sequence, we introduce a fully automatic method to detect and label those electrodes in MR images. Our method was tested on 8 subjects wearing a 64-channel EEG cap. This automated method showed an average detection accuracy of 94% and the average position error was 3.1 mm. These results suggest that the proposed method has potential for determining the position of the electrodes during simultaneous EEG/fMRI acquisition with a very light cost procedure.

Introduction

Electroencephalography (EEG) measures the electrical potential generated by the neuronal activity over the scalp with electrodes placed on the surface of the scalp [START_REF] Buzsáki | The origin of extracellular fields and currents -EEG, ECoG, LFP and spikes[END_REF], [START_REF] Murakami | Contributions of principal neocortical neurons to magnetoencephalography and electroencephalography signals[END_REF], [START_REF] Petsche | On the search for the sources of the electroencephalogram[END_REF]). Usually electrodes are placed thanks to a flexible cap and positioned according to anatomical points enabling optimal covering of brain regions regardless of the size and shape of the subject's head. Currently, when acquiring EEG and functional magnetic resonance imaging (fMRI) simultaneously, the position of the electrodes is calculated according to fiducial points (anatomical points of the skull) such as inion, nasion and vertex [START_REF] Strobel | Novelty and target processing during an auditory novelty oddball: A simultaneous event-related potential and functional magnetic resonance imaging study[END_REF]). The localisation of EEG sources in the brain depends on several parameters including the position of the electrodes on the scalp. A precise knowledge of these positions is important because inaccurate information on EEG electrodes coordinates may affect EEG inverse solution [START_REF] Khosla | Spatial mislocalization of EEG electrodes-effects on accuracy of dipole estimation[END_REF]). This knowledge is even more crucial in the case of simultaneous EEG and fMRI study, when the sessions are conducted repeatedly over a long period of time. Approximations in the positioning of the electrodes are then made in each session and will give rise to important inaccuracies in the measured evoked potential [START_REF] Wood | Interpretation of evoked potentials: A neurophysiological perspective[END_REF]). As a matter of fact, magnetic resonance (MR) images and EEG need to be registered to be able to compare activations given by fMRI and by EEG. This simultaneous acquisition allows the concordance of two different kind of information, a high temporal resolution in the order of a millisecond with EEG, and a high spatial resolution in the order of millimetre with MRI.

In this article an automated and efficient method to determine EEG electrodes positions based on a specific MR sequence is presented and evaluated. Compared to other existing approaches, the proposed method does not need additional hardware (like 3D electromagnetic digitizer devices [START_REF] Adjamian | Co-registration of magnetoencephalography with magnetic resonance imaging using bite-bar-based fiducials and surface-matching[END_REF]; [START_REF] Whalen | Validation of a method for coregistering scalp recording locations with 3D structural MR images[END_REF]), artificial electrode markers [START_REF] Sijbers | Automatic localization of EEG electrode markers within 3D MR data[END_REF]) or laser scanner [START_REF] Koessler | EEG-MRI Coregistration and Sensor Labeling Using a 3D Laser Scanner[END_REF]; [START_REF] Bardouille | Improved localization accuracy in magnetic source imaging using a 3-d laser scanner[END_REF])), which might be uncomfortable for the subject if he must stay still during acquisition [START_REF] Le | A rapid method for determining standard 10/10 electrode positions for high resolution EEG studies[END_REF]) and add time to the preparation of the patient. Semi-automated electrodes localisation methods exist [START_REF] De Munck | A semi-automatic method to determine electrode positions and labels from gel artifacts in EEG/fMRI-studies[END_REF]; [START_REF] Butler | Application of polymer sensitive MRI sequence to localization of EEG electrodes[END_REF]), which require a manual fiducial landmark identification to guide co-registration without any markers 2 2

METHODS

but these approach relies on the efficiency of the accuracy of the operator. Another automated method was recently developed and shown great results with an anatomical MR image [START_REF] Marino | Automated detection and labeling of high-density EEG electrodes from structural MR images[END_REF]), however, this method is only working with a high density cap also compatible with MRI: the GES 300 from Geodesic EEG Systems. Since this kind of cap includes plastic around electrodes and contain hydrogen protons, it can be visible on T1-w image. For seek of genericity (i.e. able to operate on all types of caps when artefacts do not appear on T1-w images), we propose to make use of a MRI sequence with radial k-space sampling named UTE for Ultra-short Echo-Time. It allows to visualise the tissues with a very short T2 and T2 * , such as cortical bone, tendons and ligaments [START_REF] Holmes | MR imaging with ultrashort TE (UTE) pulse sequences: Basic principles[END_REF], [START_REF] Keereman | MRI-based attenuation correction for PET/MRI using ultrashort echo time sequences[END_REF]). This sequence is all the more interesting in our context because it enables the visualisation of the MR compatible electrodes [START_REF] Butler | Application of polymer sensitive MRI sequence to localization of EEG electrodes[END_REF], [START_REF] Springer | Three-dimensional ultrashort echo time imaging of solid polymers on a 3-tesla wholebody MRI scanner[END_REF]) on the scalp with a capability to be performed rapidly enough to not overwhelm the whole MRI protocol.

We propose a fully automated method, which provides reliable and reproducible results for the detection and labelling of a MR compatible EEG cap into the MR space.

Methods

The retrieval of the electrodes consisted in two parts; firstly, we provided a mask that includes the volume where the electrodes are located; secondly, we performed the electrode detection inside this volume of interest (VOI). Figure 1 presents a flowchart of the method's main steps. We hypothesised that electrodes would appear as spheres inside the UTE volume and it allows us to perform a Hough transform in a consistent manner across subjects.

Scalp segmentation

Several reliable scalp segmentation methods exist for T1-w imaging. Because UTE images are noisier, we performed the scalp segmentation on the T1-w images and co-registered the UTE images with the T1-w images to apply the mask. The T1-w is first registered on the UTE and the anatomical T1 image is then segmented using FSL, an open-library of analysis tools for MRI and its function BET (Brain Extraction Tool) [START_REF] References | Fast robust automated brain extraction[END_REF]; [START_REF] Popescu | Optimizing parameter choice for FSL-Brain Extraction Tool (BET) on 3D T1 images in multiple sclerosis[END_REF]). A mask of the scalp is computed from the segmentation. Since electrodes are located around the head of the subject, the scalp mask is dilated toward the periphery in order to isolate this layer. What is outside the dilated mask is subtracted in order to isolate only the layer where the electrodes are located.

Detection of electrodes with the Spherical Hough transform

A 3D Hough transform was used to segment the electrodes inside the VOI. Hough transform is typically used to detect circles or lines in 2-dimensional data sets, but was recently extended to detect spheres in 3-dimensional data sets [START_REF] Xie | Magnetic resonance histology of age-related nephropathy in the Sprague Dawley rat[END_REF]; [START_REF] Borrmann | The 3d hough transform for plane detection in point clouds: A review and a new accumulator design[END_REF]). As the shape of an electrode can be assimilated to a sphere, the Spherical Hough Transformation algorithm seemed particularly well adapted to this task. The VOI image is first smoothed using a Gaussian kernel, with a FWHM (Full Width at Half Maximum) adapted to the size of the electrode (10 mm) in order to reduce the noise of the image while saving electrode information. Then, the Hough algorithm is performed and provides a list of n potential electrodes, D = [d 1 , . . . , d n ]. Figure 2 shows an example of such detections on a 2D slice of the VOI. Because the VOI includes also anatomical structures (nose, ears) and noise (artefacts due to the cap or gel), the number of potentially detected electrodes is substantially higher than the number of "true" electrodes N , in our case 64.

4

Selection of detected electrodes

The detected electrodes are then filtered to get rid of the potential false detections given by the Hough transform. A 64 electrodes spherical EEG template p j (1 ≤ j ≤ 64) ∈ P was given by the cap manufacturer, indicating theoretical positions of every electrodes relatively to each other. Due to the non-sphericity of the head and the elastic deformations of the cap, these positions are not sufficient enough to give a reliable detection by itself. However, this template will be used to identify outliers in our detections. This spherical template is registered onto the detected electrodes from previous section, through the Iterative Closest Point (ICP) algorithm, a well-known algorithm for registering two-cloud of points [START_REF] Besl | A method for registration of 3-d shapes[END_REF]; [START_REF] Chen | Object modelling by registration of multiple range images[END_REF]). The algorithm takes a first point cloud which will be kept fixed, while the other one will be spatially transformed in order to best align the reference. The goal is to iteratively minimise a metric error, usually the distance between the two sets of points, by modifying the transformation applied to the source.

In our case, the ICP will find the optimal rotation, translation and scale to fit the data point set D obtained with the Hough transform and the model point P . The algorithm is divided into 2 steps. The first step consists in estimating correspondences between the two set of points. During this step, for each point p j , in the reference set P , the closest point d i of the detected points set D is computed. This point will be noted c j and therefore defined as follows:

c j = arg min d 1 ,...,dn dist(d i , p j ), ∀j ∈ [1, . . . , N ].
(1)

The second step consists in computing the similarity transform that best aligns every c j to the corresponding p j . The minimisation is expressed by:

(R * , S * , t * ) = arg min R,S,t j∈[1,...,N ] c j -SR p j -t 2 , ( 2
)
where R is a rotation matrix (3 × 3), t is a translation vector (3 × 1) and S is a scale matrix (S = s * Id, 3 × 3). The ICP runs until convergence. The registered template P can then be written 5 5

METHODS

as:

p j = SR p j + t. (3)
Once the ICP is completed, a two-part filtering phase is implemented. The first one consists in taking the closest point of the Hough transform data set; for each of the N electrodes of the registered model P , the closest detected point c j is selected. Unselected points are discarded and, after this first filtering step, the number of electrodes is therefore equal to N , the total number of electrodes desired (64 in our case). Figure 3 illustrates the impact of this step. For the second and final step, all points c j , which are too far from the closest point of the template P , are removed. A threshold equals to four times the Median Absolute Deviation (MAD) of all distances is applied. For each removed point, a replacement is determined by a new detection from the local maxima on the VOI image around the theoretical position given by the registered template (cf. Figure 4). The new data set D is obtained and the N electrodes are then labeled using the template.

6 6

Validation of the method

A manual selection of the electrodes positions was done on the UTE sequence and the quality of our detection was assessed using this manual selection as a ground truth. Instead of selecting the center of each electrode in a 3D image, we choose to use a more convenient procedure for the manual detection. Following [START_REF] Butler | Application of polymer sensitive MRI sequence to localization of EEG electrodes[END_REF], the manual detection was performed by picking up the Cartesian position (x i , y i , z i ) of each 64 electrodes for each subject on a pancake view, which is roughly a 2D projection of the scalp (de Munck et al. [2012]).

The performance indicators of our automated detection will be the position error (PE) and the positive predictive value (PPV). The position error is the average Euclidean distance between each pair of electrodes (the manually selected one, considered as the ground truth, and the detected one) and the PPV is the percentage of electrodes that have been well detected. We considered that a detected electrodes is well localised when the PE is below 10 mm, which corresponds to the diameter of the electrode [START_REF] Kavanagk | Evaluation of methods for threedimensional localization of electrical sources in the human brain[END_REF]).

We also compared the performance of our method against a more traditional semi-automatic one: five fiducial points were selected manually and the spherical template was adjusted to these points [START_REF] Towle | The spatial location of EEG electrodes: locating the best-fitting sphere relative to cortical anatomy[END_REF]). This method, although not recent, is still used by many studies (e.g. [START_REF] Thornton | Sensorimotor activity measured via oscillations of EEG mu rhythms in speech and non-speech discrimination tasks with and without segmentation demands[END_REF][START_REF] References | Trait related sensorimotor deficits in people who stutter: An EEG investigation of µ rhythm dynamics during spontaneous fluency[END_REF], [START_REF] Ge | Temporal-Spatial Features of Intention Understanding Based on EEG-fNIRS Bimodal Measurement[END_REF]).

Subjects and EEG equipment

After IRB approval, eight healthy volunteers provided written informed consent to take part in the study. They all underwent a simultaneous EEG/fMRI examination (fully described in [START_REF] Mano | How to build a hybrid neurofeedback platform combining EEG and fMRI[END_REF]). EEG was acquired using two 32-channel MR compatible amplifiers (actiCHamp, Brainproduct, Gilching, Germany) and a cap providing 64 Ag/AgCl electrodes positioned according to the extended 10-20 system and one additional ground electrode. Electrodes are attached to small cups with inner diameter of 10 mm and 4 mm height, inserted in the cap and filled with gel to minimise the contact impedance. All subject wore a large (circumference between 56-58 cm) MR compatible cap from Brainproduct (Gilching, Germany) and a particular attention was given to its positioning according to standard fiducial points.

UTE sequences parameters

All MR data were collected on a 3T Siemens Verio MR scanner (VB17, Siemens Healthineers, Erlangen, Germany). Specifically, the UTE sequence using 3D radial k-space sampling was performed with the following parameters: repetition time (TR) = 3.45 ms, echo time (TE) = 0.07 ms, flip angle (FA) = 14 • and voxel size 1.33 × 1.33 × 1.33 mm 3 . A 3D T1 MPRAGE was also performed: TR = 1900 ms, TI = 900 ms, TE = 2.26 ms, FA = 9 • and voxel size 1 × 1 × 1 mm 3 . Two additional UTE sequences with lower sampling resolution were acquired in order to decrease the acquisition time and to investigate the impact on electrodes detection. To reduce the acquisition time, the number of spokes has to decrease; from 60000 spokes (60K) for the original, to 30000 (30K) and 15000 (15K) spokes for the additional ones. The UTE acquisition time goes down from 5 min 35 s to 2 min 47 s and 1 min 23 s. A comparison between these acquisitions is shown in Figure 5.

Results

The creation of an image (VOI) containing only the information related to the electrode allows to remove external noise while protecting the information related to the electrode. This image enables robust detection of the position of the electrodes for all subjects. Furthermore, since our method always detects exactly N (64 in our case) electrodes, the number of false negatives (missed electrodes) will automatically be equal to the number of false positives (wrongly detected electrodes). Table 1 presents the mean position error (PE), the standard deviation of the PE and the maximum PE of our detections for each of the eight subjects. The max PE reflected a high difficulty to detect the electrodes near anatomical parts or in posterior regions where the head apply a pressure on the EEG cap inside the MRI. Our UTE-based electrode detection showed an average PE of 3.1 mm for all subjects. The detection accuracy, represented by the positive predictive value (PPV), is also shown and corresponds to the percentage of electrodes correctly found. The average PPV for all subjects was 94.22%. The PPV is the percentage of electrodes that have been detected. We consider that an electrode is well localised when the PE is below 10 mm, which represents the diameter of an electrode. The mean PE on all subject is equal to 3.1 mm and the mean PPV to 94.22%.

We then compared the performance of our method with the semi-automatic one presented in section 3.4 (FID). The PE and PPV were calculated in the same way. The results are shown in Table 2 and Figure 6 shows, for each subject, a comparison of the PEs obtained by the two methods. The mean PE on all subject is equal to 7.7 mm and the mean PPV to 79.41%. Moreover, for every subjects, our method produced smaller PE and better PPV. A paired t-test was computed between the two PEs sets and a significant difference was obtained (p<0.0001). 9 9 ) for each subject (S1-S8) for semi-automated electrodes detection based on manual delineation of fiducial landmark (FID). The PPV is the percentage of electrodes that have been detected. We consider that an electrode is well localised when the PE is below 10 mm, which represents the diameter of an electrode. The mean PE on all subject is equal to 7.7 mm and the mean PPV to 79.41%. Finally, we investigated the impact of lower sample UTE sequences, which allow reducing the acquisition time, on electrode detection. We tested two others UTE sequence (cf. Section 4.2). We applied our detection method on the three different UTE images and compared the quality of the detections. Table 3 reports the mean PE and mean PPV obtained for the three UTE sequences on seven subjects (the first subject did not receive the additional sequences). As expected, the mislocalisation, as well as the position error, increase according to the decrease of the sampling. However, our results are still clearly better than the semi-automatic one for the 30k sequence (half the acquisition time than the original one) and are slightly better for the fastest sequence.

RESULTS

Subjects

Discussion

We have proposed an automated method for detecting and labelling EEG electrodes based on UTE MR images without using any external sensors. Previous results indicate that a localisation technique using electromagnetic digitisation technology is time-consuming [START_REF] Dalal | Consequences of EEG electrode position error on ultimate beamformer source reconstruction performance[END_REF]) and others techniques such as 3D digitisation can be affected by errors of registration and projection of EEG electrodes on the head model. We have shown that our method offers constant and precise results. Moreover, the proposed method provides the position of the electrodes directly into the MR-space, which is crucial in case of simultaneous EEG/fMRI acquisitions.

Furthermore, for seek of genericity, the proposed method is able to operate on all types of caps and does not need specific electrodes, unlike a recent work from [START_REF] Marino | Automated detection and labeling of high-density EEG electrodes from structural MR images[END_REF] for example. To the best of our knowledge, this is first automated electrodes detection method implying non-visible electrodes on anatomical MR sequence.

The method presented here requires only an additional sequence (the UTE acquisition sequence) in the experimental protocol. This acquisition takes from 1 to 5 min. From our experiments, a good compromise between acquisition time and detection quality can be achieved with a 2 or 3 min sequences. Further optimisation of the sequence parameters could enable an improvement of the images without increasing the acquisition time.

Conclusion

We presented a method to automatically detect and label EEG electrodes during an EEG/fMRI acquisition. We used a UTE MR sequence to obtain electrodes positions on a MR-volume. This method only has for additional cost the acquisition time of the UTE sequence in the MR protocol. We have demonstrated that our method achieves a significantly more accurate electrode detection compared to a semi-automatic detection one that is more commonly used during EEG/fMRI protocols. For future research, since the proposed method can be totally automated and does not require complex processing, this technique may be used to extract the position and the label of the electrodes in real time. Indeed, this technique is interesting for applications requiring immediate knowledge of the position of the electrodes. We believe this method will be useful to improve the fusion of EEG and fMRI signals.

11

MRI data acquisition was partly supported by the Neurinfo MRI research facility from the University of Rennes I. Neurinfo is granted by the the European Union (FEDER), the French State, the Brittany Council, Rennes Metropole, Inria, Inserm and the University Hospital of Rennes. The authors would like to thank their collaborators at Siemens Healthineers, in particular Vladimir Jellus and P. Speier for their valuable contributions to the development and testing of the prototype Ultra Short TE for this project.

Figure 1 :

 1 Figure 1: Steps for the extraction of the Volume Of Interest (VOI). An outskin mask is performed from the T1 image (1), then a dilation and a removal of the mask is performed (2) in order to obtain the layer where the electrodes are located. Finally, the UTE image is masked by the dilated mask (2) which gives us the VOI (3).

Figure 2 :

 2 Figure 2: Example of Hough transform detection (red dots) on the VOI smoothed image. Hough transform detects also anatomical parts (arrow), which will be excluded in the filtering steps (cf. section 3.3).

Figure 3 :

 3 Figure 3: Example of outliers removal in potential electrodes data set D with the ICP algorithm. The dataset D is represented in red on the left along with the registered template P in purple. The data set obtained after the first filtering step is in red on the right. Outliers are mostly due to external anatomical parts or noise not taken in account during the segmentation. These outliers are discarded by the filtering step because they are too far from P .

Figure 4 :

 4 Figure 4: Cross section of the VOI image. Green points are corresponding to the template data set P , blue points to the maximum local detection and the red one are the outliers from D. The second and final filtering step consists in replacing any point from the Hough data set too far from the registered template P . The substituted point comes from a detection by local maxima, closest to the template P .

Figure 5 :

 5 Figure 5: Example of UTE images with different sampling. The image quality as well as the acquisition time decrease linearly according to the sampling. Acquisition time for 1) 5 min 35 s, 2) 2 min 47 s, 3) 1 min 23 s.

Figure 6 :

 6 Figure 6: Position Error (PE) for UTE-based electrodes detection method (UTE) and the semi-automatic method based on fiducial points (FID). Box-plots for the eight subjects are shown.

Table 1 :

 1 Position error (PE) and positive predictive value (PPV) for each subject (S1-S8) for UTE-MR electrodes detection.

	Subjects Mean PE (mm) Std PE (mm) Max PE (mm) PPV (%)
	UTE S1	2.73	2.83	17.32	95.38
	S2	2.41	2.38	12.43	96.92
	S3	2.66	2.43	13.62	95.38
	S4	4.10	3.97	21.3	89.23
	S5	3.32	3.38	14.45	90.7
	S6	2.85	3.24	17.57	95.38
	S7	2.76	2.51	15.90	96.92
	S8	3.69	5.13	26.24	93.84

Table 2 :

 2 Positive predictive value (PPV) and position error (PE

		Mean PE (mm) Std PE (mm) Max PE (mm) PPV (%)
	FID S1	7.60	3.09	15.31	76.92
	S2	7.79	2.85	14.63	72.30
	S3	6.15	2.71	13.45	90.76
	S4	7.70	3.41	17.78	73.84
	S5	6.08	2.78	16.36	90.76
	S6	6.54	3.40	17.47	87.69
	S7	6.88	3.39	17.10	87.69
	S8	12.55	6.36	30.6	55.38

Table 3 :

 3 Mean of position error (PE) and mean positive predictive value (PPV) for three different sampling resolutions of the UTE sequence. Shorter acquisition time implied lower SNR and lower detection accuracy. Results are still better than the semi-automatic method.