Alexis Piedrafita

Stéphane Balayssac

Audrey Casemayou

Jean-Sébastien Saulnier-Blache

Alexandre Lucas

Jason Iacovoni

Benjamin Breuil

Dominique Chauveau

Stéphane Decramer

Myriam Malet-Martino

Aurélie Monié

Sophie Franceschi

Myriam Malet- Martino

Mathieu Delample

Emile Perez

Jean-Christophe Garrigues
email: garrigues@chimie.ups-tlse.fr

Study of rapeseed oil gelation induced by commercial monoglycerides using a chemometric approach

Keywords: Oleogel, Acylglycerides, Chemometrics, QSPR, Artificial neural network, GC, NMR, Rapeseed oil

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

Introduction

Fats and oils are largely used as main ingredients in food preparations. Because of their high composition in saturated fatty acids (FAs), the fats are solids at room temperature and this useful property is responsible for the sensory attributes of numerous products [START_REF] Marangoni | Structure and functionality of edible fats[END_REF]. Due to the development of cardiovascular diseases as well as environmental concerns (in the case of palm oil), hard fats are increasingly disliked by consumers and manufacturers [START_REF] Hooper | Reduction in saturated fat intake for cardiovascular disease[END_REF][START_REF] Khatun | Fatty acid composition, fat deposition, lipogenic gene expression and performance of broiler fed diet supplemented with different sources of oil[END_REF]. Besides, some saturated FAs have been implicated in adverse effects to health, such as trans-FAs [START_REF] Booker | Trans fatty acids and cardiovascular health: Translation of the evidence base[END_REF]. For all of these reasons, there is nowadays a real demand for solid lipids and consequently there is growing interest in oil structuring alternatives for food industry [START_REF] Doki | Virtual bioequivalence for achlorhydric subjects: The use of PBPK modelling to assess the formulationdependent effect of achlorhydria[END_REF][START_REF] Martins | Edible oleogels: An opportunity for fat replacement in foods[END_REF][START_REF] Martins | Oleogels for development of health-promoting food products[END_REF][START_REF] Patel | CHAPTER 1 Oil Structuring: Concepts{,} Overview and Future Perspectives[END_REF][START_REF] Patel | Applications of fat mimetics for the replacement of saturated and hydrogenated fat in food products[END_REF][START_REF] Pehlivanoğlu | Oleogels, a promising structured oil for decreasing saturated fatty acid concentrations: Production and food-based applications[END_REF][START_REF] Puscas | Oleogels in food: A review of current and potential applications[END_REF].

Oleogelation is a promising technology used in food chemistry to design soft matter structures with new properties in which the liquid oil is entrapped through the formation of a threedimensional network. The oleogelators used to form oleogels can be classified into different structuring networks: self-assembly fibers (12-hydroxystearic acid, lecithins with βsitosterol), crystalline particles (waxes, monoglycerides (MGs)), and polymeric gelators like cellulose derivatives such as ethyl cellulose, carboxymethyl cellulose and hydroxyl propyl methyl cellulose [START_REF] Doan | Internal and external factors affecting the crystallization, gelation and applicability of wax-based oleogels in food industry[END_REF][START_REF] Han | Rheological properties of organogels developed by sitosterol and lecithin[END_REF][START_REF] Meng | Oleogels from sodium stearoyl lactylatebased lamellar crystals: Structural characterization and bread application[END_REF][START_REF] Naeli | Development of innovative ethyl cellulose-hydroxypropyl methylcellulose biopolymer oleogels as low saturation fat replacers: Physical, rheological and microstructural characteristics[END_REF][START_REF] Palla | Multi-objective optimization of high oleic sunflower oil and monoglycerides oleogels: Searching for rheological and textural properties similar to margarine[END_REF][START_REF] Rogers | Nanostructuring fiber morphology and solvent inclusions in 12-hydroxystearic acid / canola oil organogels[END_REF]. Among the food grade oleogelators, MGs are very interesting for their gelation properties, but also from an economic point of view (Co & Marangoni, 2012;[START_REF] Yu | Development of a food-grade organogel with high bioaccessibility and loading of curcuminoids[END_REF]. These gelators also have the ability to structure both oily and aqueous phases [START_REF] Sein | Rheological characterization, crystallization, and gelation behavior of monoglyceride gels[END_REF]. For all these reasons, oil gelation by MGs for food applications has been the subject of numerous studies [START_REF] Giacomozzi | Muffins Elaborated with Optimized Monoglycerides Oleogels: From Solid Fat Replacer Obtention to Product Quality Evaluation[END_REF][START_REF] Kouzounis | Partial replacement of animal fat by oleogels structured with monoglycerides and phytosterols in frankfurter sausages[END_REF][START_REF] Palla | Multi-objective optimization of high oleic sunflower oil and monoglycerides oleogels: Searching for rheological and textural properties similar to margarine[END_REF]. Below a critical gelation temperature, MGs can form compact reverse lamellar phases which will lead to lamellae or platelets. The latter will then assemble into crystalline clusters which will be able to percolate and form a three-dimensional network that will entrap the oil (gelation) [START_REF] Marangoni | Edible oleogels : structure and health implications[END_REF]. Gelation is therefore based on a crystallization phenomenon, and anything that can affect it (ordering or disordering) will affect the gel and its strength. This phenomenon also depends on the nature of the oil [START_REF] Ferro | Role of the oil on glyceryl monostearate based oleogels[END_REF], the MGs concentration [START_REF] Ojijo | Effects of monoglyceride content, cooling rate and shear on the rheological properties of olive oil/monoglyceride gel networks[END_REF], and their intrinsic structure [START_REF] Cerqueira | Structural and mechanical properties of organogels: Role of oil and gelator molecular structure[END_REF]. The long chain length and the degree of unsaturation play an important role in the gelation process [START_REF] Vereecken | Comparing the crystallization and polymorphic behaviour of saturated and unsaturated monoglycerides[END_REF].

Commercial MGs for industrial use (even relatively pure) are synthesized from natural products (vegetable oils) and remain complex mixtures composed of different chain lengths and unsaturations that are difficult to purify [START_REF] Nitbani | Preparation of fatty acid and monoglyceride from vegetable oil[END_REF]. In addition, they are generally mixed with variable amounts of residual diglycerides (DGs) (and their respective isomers), triglycerides (TGs) and free FAs; all compounds present different chain lengths and degrees of unsaturation. Consequently, all these ingredients and their structures can mutually influence the gelation process. The mechanisms involved in oil gelation by DGs oleogelators involve several steps of crystallization, which are influenced by their isomer structures. In specific mixtures containing two isomers of 1,2-and 1,3-DGs, significant differences in the strength of the hydrogen bonds of the hydroxyl groups and the aliphatic chains, related to the specific three-dimensional arrangements, are observed: V-shaped FA chain arrangement for 1,3-DGs and hairpin-shaped conformation for 1,2-DGs [START_REF] Lo | Diacylglycerol oilproperties, processes and products: A review[END_REF]. MG structures will also influence crystallization mechanisms based on their saturated and unsaturated FA composition and their chain lengths [START_REF] Da Silva | Crystallisation of monoacylglycerols and triacylglycerols at different proportions: Kinetics and structure[END_REF]. Different pure MGs, and commercial food-grade mixtures of MGs, have also been studied for their gelation properties with vegetable oils, showing a link between the mixture's composition and rheological properties [START_REF] López-Martínez | Comparing the crystallization and rheological behavior of organogels developed by pure and commercial monoglycerides in vegetable oil[END_REF]. To date, studies concerning the gelation of oil with commercial MG products have not been fully described considering the positive or negative interactions related to the structure and the composition of all the species present in the mixture: TGs, MGs, 1,2-and 1,3-DGs, and FA composition. The mechanisms involved in structuring acylglyceride (AG) oleogelators within the oil matrix are complex and not fully understood. The first step is the solubilization of the oleogelator molecules by the liquid oil phase at high temperature, followed by the crystallization of the oleogelators while cooling, reducing the mobility of the oil within a structured network. The properties of this structured network can be studied using rheological measurements [START_REF] Patel | Rheological Profiling of Organogels Prepared at Critical Gelling Concentrations of Natural Waxes in a Triacylglycerol Solvent[END_REF]. The oleogelator aggregates can act as building blocks to trap the liquid oil into a gel or a gel-like structure, depending on the nature and concentration of the oleogelators. The structure of the building blocks is stabilized by weak interactions such as H-bonding or dipole. Depending on the nature of these interactions, different types of aggregates will be formed at the microscopic scale, resulting in strong or weak structures, and producing viscous liquids or gels at the macroscopic level that can be studied in the linear viscoelastic region [START_REF] Fayaz | Comparative Study on the Ability of Different Oleogelators to Structure Sunflower Oil[END_REF]. The analysis of complex mixtures of AGs is a difficult challenge, and the analytical strategy for characterizing the composition of mixtures must be defined between targeted and nontargeted methods. Recently, we have developed a targeted method based on high-resolution mass spectrometry and liquid chromatography to identify each component of acylpolyglycerol, transposable to commercial MGs characterization [START_REF] Garrigues | Analysis of complex mixtures of polyglycerol fatty esters using liquid chromatography and high-resolution mass spectrometry: Retention, structural and relative composition study[END_REF]. Due to the unknown origin of the commercial MGs, we have chosen in the present study to evaluate a non-targeted analytical methodology based on proton Nuclear Magnetic Resonance

( 1 H NMR) and gas chromatography (GC). The composition of the commercial MGs was determined by quantification of AGs families, determination of the different FAs present in these mixtures and the calculation of the ratio between these different analytical descriptors. This analytical methodology was then used to study the gelation of rapeseed oil with different concentrations of commercial MGs. This oil was specifically chosen due to its health and nutritional benefits, with its high concentration in both oleic and alpha-linolenic acids [START_REF] Vingering | Fatty acid composition of commercial vegetable oils from the French market analysed using a long highly polar column[END_REF][START_REF] Xie | Profiling and quantification of lipids in cold-pressed rapeseed oils based on direct infusion electrospray ionization tandem mass spectrometry[END_REF]. The classification of the commercial MGs was achieved through the use of a quantitative structure/composition properties relationship (QSPR) by artificial neural networks (ANNs) and multivariable linear regressions (MLRs) supported by analytical parameters obtained by GC and 1 H NMR. Several ratios between specific chemical classes of lipids were introduced, in addition to analytical parameters such as structure and composition descriptors. All these factors were used to relate the gelation of rapeseed oil to molecular data, evaluated by rheological measurements (G'). For a wider investigation of the gelation properties, we used thirteen different commercial MGs mixtures.

Materials and methods

Chemicals and reagents

Rapeseed oil (density at 20 °C = 0.92; viscosity = 80 cP) was used as the liquid oil for all dispersions with commercial mixtures of mono and diglycerides (MDGs) and was purchased in a convenience store (METRO ® ). The commercial MDGs used in this study come from several providers: MDG 0096, MDG 5616, MDG 0298 and MDG 0051 (Palsgaard), Lamemul K 2000 K (BASF), Multec Mono 9202 (Puratos) and the medium-chain MDGs (type I & type II), glycerol dibehenate, glycerol monostearate as well as hemi-synthetic glycerides (HSglycerides) (Stearinerie Dubois). Type I medium-chain MDGs are richer in C10:0 FA than Type II.T These distilled MDGs are obtained from edible vegetable oil. All these emulsifiers are specifically purified and produced for use in food applications. Various reactants were used for the specific analyses. Sodium methoxide, sodium hydroxide, boron trifluoride, isooctane and methanol were provided by Sigma Aldrich. The solvent for 1 H NMR analyses, deuterated chloroform, CDCl3 (>99.5%), was purchased from Eurisotop.

Analysis of the commercial AGs mixtures by GC

All MDG characterizations were conducted at ITERG (Institut des Corps Gras, Pessac, France). All analyses were performed in triplicate.

Procedure for medium-chain MDGs

Transesterification was performed with sodium methoxide (5% w/w). Methyl esters were analyzed on a GC apparatus (Agilent 6800) equipped with a split/splitless injector, an FID detector and a CP-Sil 88 polar column (100 m/0.25 mm/0.2 μm). Injector temperature: 250 °C; injection ratio: 1/10; FID temperature: 275 °C. The temperature program was: 60 °C for 5 min; temperature ramp 1: 15 °C/min until 165 °C; temperature ramp 2: 2 °C/min until 225 °C followed by a constant temperature for 20 min as the final step. The results were obtained after applying a specific response factor to each FA methyl ester. Quantification was carried out using the C11 methyl ester as the internal standard. All values are expressed with a minimum and a maximum deviation of 0.5 and 3.5%, respectively.

Procedure for long-chain MDGs

Fat saponification was achieved by a reflux with 0.5M NaOH in methanol and esterification was performed with 14% boron trifluoride (BF3) in methanol. The samples were then dissolved in isooctane and analyzed on a CP-Sil 88 polar column (100 m/0.25 mm/0.2 μm). Samples were injected according to the following program: injection pressure: 225 KPa; injector temperature: 250 °C, detector: 275 °C, injection ratio: 1/100. The temperature program was: oven at 60 °C for 5 min, temperature ramp 1 at 15 °C/min until 165 °C, temperature ramp 2 at 2 °C/min for 1 min and a constant temperature at 225 °C for 20 min as a final step. The results were obtained after applying a specific response factor to each FA methyl ester. Quantification was carried out using the C17 methyl ester as the internal standard. All values are expressed with a minimum and a maximum deviation of 0.5 and 3.5%, respectively.

Analysis for the commercial acylglyceride mixtures by 1 H NMR

A known amount of each MDGs was diluted with 700 μL of CDCl3, mixed and then transferred to a 5-mm NMR tube. The 1D 1 H NMR experiment was acquired on a Bruker Avance 500 spectrometer (Bruker Biospin) equipped with a dual 1 H-13 C TCI cryoprobe at 25 °C. A conventional NMR pulse sequence with a relaxation delay of 6.54 s, a flip angle of 30° and an acquisition time of 5.46 s was used. NMR data were processed using the Topspin 4.08 software (Bruker Biospin) with one level of zero-filling and Fourier transformation after multiplying FIDs by an exponential line-broadening function of 0.3 Hz. Phase adjustment and polynomial baseline correction were done manually on each spectrum. All experiments were performed in duplicate and the deviation between experiments never exceeded 5%.

Integration was performed on the CH group of the glycerol moiety, at 3.9 ppm (1-MGs), 4.9 ppm (2-MGs), 5.1 ppm (1,2-DGs), 4.1 ppm (1,3-DGs) and5.3 ppm (TGs). The percentage of each type of molecule was estimated by dividing the area of its NMR signal by the total area of all integrations.

Preparation of the gels

Mixtures of MDGs at different concentrations (from 2.5 to 20 g/100 g) were prepared by dispersing MDGs (in a liquid, paste or powder form) in rapeseed oil. All the mixtures were warmed at 80 °C for 1 hour. At this temperature, and after manual mixing, all solutions appeared isotropic. The samples were stored at 20 °C and some of them became more or less turbid upon cooling, obtaining numerous textures from liquid to solid. Samples were then analyzed by rheology after a crystallization equilibrium period (i.e. 1 week).

Characterization of the mixtures by rheological measurements

Oscillatory experiments were performed at 20 °C with an ARG2 rheometer (TA Instruments) on MDG/oil mixtures. All experiments were conducted in a common linear regime valid for high and low textures, as previously described for the study of hard and soft oil gels [START_REF] Fayaz | Comparative Study on the Ability of Different Oleogelators to Structure Sunflower Oil[END_REF] at a deformation amplitude of 0.1%, and at a frequency of 0.5 Hz.

Parallel-plate geometry was used with a gap of 1 mm and rough surfaces in order to avoid wall slipping. The MDG/oil mixtures were placed in the cell of the rheometer very carefully; sample compression during loading was the minimum possible and was always done at the smallest rate. Excess sample was gently removed from the edge of geometry and the whole apparatus was covered during the analysis (global time: 2.5 min with a storage modulus acquisition (G'; Pa) every 5 s). Temperature within the sample was controlled with a precision of 0.1 °C via a thermoelectric Peltier module. All rheological measurements were conducted in triplicate and the deviation between experiments never exceeded 5%.

Analytical database

An analytical database was constructed for each commercial AG mixtures by integrating the GC results: FA levels: C6:0; C8:0; C10:0; C12:0; C14:0; C15:0; C16:0; C16:1; C17:0; C17:1; C18:0; C18:1; C18:2; C18:3; C20:0; C20:1; C22:0; C22:1; C24:0; saturated FAs; monounsaturated FAs; polyunsaturated FAs and trans FAs. From these raw GC data, specific analytical parameters were calculated: intermediate FAs ( 6 long-chain/very long-chain FAs. Specific ratios were also calculated from the raw 1 H NMR data, and to take into account the low concentrations of certain AGs classes, the inverse ratios were calculated: MGs/DGs; MGs/TGs; DGs/MGs; DGs/TGs; TGs/MGs; TGs/DGs; 1-MGs/2-MGs; 2-MGs/1-MGs; 1-MGs/DGs; 1-MGs/TGs; 2-MGs/DGs; 2-MGs/TGs; 1,2-DGs/1,3-DGs; 1,3-DGs/1,2-DGs; 1,2-DGs/MGs; 1,3-DGs/MGs; 1,2-DGs/TGs; 1,3-DGs/TGs; 1-MGs/1,2-DGs; 1-MGs/1,3-DGs; 2-MGs/1,2-DGs; 2-MGs/1,3-DGs. If the denominator of one of these ratios is zero, the ratio was not calculated and a null value was assigned to it.

Chemometrics and statistical analysis

All data obtained in this study were normalized by dividing each variable by the greatest value of the column. The QSPR study by ANNs was carried out in MATLAB R2020b (MathWorks) with the deep learning toolbox V.14.1. An ANN was constructed with 64 input neurons, connected to each of the 64 analytical parameters and one output neuron connected to the G' value measured for the 5 concentrations of each of the 13 commercial AGs.

Feedforward backpropagation was used as a learning algorithm, with a sigmoid transfer function. No hidden layer was used to extract, after 25 learning cycles, the weight attributed by the ANNs to each of the 64 analytical parameters.

Multiple linear regression models were built in MATLAB R2020b, with the statistics and machine learning toolbox V.12.0, using the key analytical parameters identified by the ANNs. For each model, the coefficient of determination (R²), the root-mean-square error (RMS) and the p-value from Fisher's test were calculated from the normalized values.

Principal component analysis (PCA) works by using combinations of variables to describe the variance of a matrix by projecting it on different components. Eigenvectors represent each variable in a correlation matrix for all components, and the commercial AGs can be plotted on the corresponding score plot. PCAs were performed in MATLAB R2020b, with the statistics and machine learning toolbox V.12.0, using the key analytical parameters identified by the ANNs as variables for each concentration in commercial AGs mixtures to construct the components. From all the components of the correlation matrices, 3 were included in the analyses to perform classification of the 13 commercial AGs mixtures and 3D plots were defined to represent the correlation score for all concentrations used in the gelation process.

Results and discussion

Analytical database

First of all, the 13 commercial AGs mixtures composed of 11 known oleogelators and 2 negative controls (C8-C10 where no gelation occurs) were analyzed by GC and 1 H NMR to identify and quantify each lipid classes. With this information, a structural database was built with 64 analytical parameters (Table S1, supplementary data). To construct Table S1, we selected the analytical parameters suitable for describing the structure of AGs mixtures, grouped into 4 families: the lipid class of FAs, the composition of FAs, the amount of TGs and the different isomers of DGs and MGs and the ratios calculated between these different structural descriptors. The ratios calculated allow a relative comparison of the commercial AGs mixtures and any unknown minor impurities are not taken into account for the QSPR study, which may, in this case, ignore them (Cherkasov et al., 2015).

The composition of the FAs classes monounsaturated, long, intermediate and the corresponding AGs varies widely depending on the mixtures investigated. Some analytical descriptors (Table S1, supplementary data) show high variability in the composition of the oleogelators,with saturated FAs varying from 7.9 to 99.9 g/100g, mono-unsaturated FAs ranging from 0.06 to 83.3 g/100g, intermediate FAs (C6-C10) varying from 0.1 to 99.6 g/100g and long-chain FAs (C12-C22) varying from 0.0 07 to 99.8 g/100g. This composition is responsible for a variable gelation effect on rapeseed oil when varying the concentration in the media, as shown in Table 1. A specific codification was applied in this case: 0 means no gelation; 0.5: liquid-solid equilibrium with flowing and 1: strong gelation with a solid appearance and no oil clearance after crystallization.

Macroscopic analyses of the inverted tube are obtained after 24 hours by visually observing the structure of the oil formulated with increasing concentrations of the 13 oleogelators. 1). Soft gels were defined and studied by Fayaz et al. [START_REF] Fayaz | Comparative Study on the Ability of Different Oleogelators to Structure Sunflower Oil[END_REF] with certain types of MGs for which macroscopic observations showed the presence of a residual liquid oil phase in the structured sunflower oil.

This observation highlights the complexity of the systems and the importance of the mixing effects, which are still difficult to study.

QSPR by ANNs

To study the contribution to gelation of all the identified parameters, we used the analytical database (Table S1, supplementary data) as input descriptors of the AGs mixtures in a QSPR study. The first step was done by analyzing the weights assigned by ANNs to each analytical parameter in relation to the storage modulus (G') of the rapeseed oil media formulated with AGs mixtures and measured by stress-controlled rheology. Table 2 presents the results of the G' measurements for the 13 AGs mixtures studied, based on their concentration in rapeseed oil. Table S2 (supplementary data) gives the key parameters with weight above 0.6 for all the concentrations of AGs mixtures in rapeseed oil and Fig. S1 (Supplementary data) summarizes all the results with positive and negative contributions to G' according to the structural data measured by analytical parameters. For concentrations between 5-15 wt.%, The C16:0 FA chain length of the AGs showed high weights for all the concentrations studied (2.9, 5%; 2.5, 5%; 5.6, 10%; 3.7, 15%; Table S2, supplementary data).

AGs with C18:0 FA showed a negative contribution at low concentrations of AGs mixtures (-1.5, 2.5%; -1.0, 5%; Table S2, supplementary data), no contribution for 10 and 15 wt.% and a positive contribution for the highest concentration tested (1.5, 20%; Table S2, supplementary data). This observation is in agreement with the fact that crystallization and melting temperature of MGs increase with chain length [START_REF] Hagemann | Effects of chain length, conformation and α-form packing arrangement on theoretical monoacid triglyceride β′-forms[END_REF]. Nevertheless, it is more difficult to understand the impact of two carbons of difference (between C16 and C18) on the gel strength as a function of the concentration. This difference could be due to a more complex crystallization polymorphism, however, as described for MGs with a chain length of 18 carbon atoms and longer [START_REF] Hagemann | Effects of chain length, conformation and α-form packing arrangement on theoretical monoacid triglyceride β′-forms[END_REF]. In order to ensure a classification of the gelation properties of these compounds, it is important to correlate these FAs to the AGs classes (MG, DG, TG) and to their positions. It appears from this study (Fig. S1 and Table S2, supplementary data) that the 1,3-DGs/1,2-DGs ratios have a positive effect (from 0.6 to 2.8; Table S2, supplementary data) for concentrations at 2.5%, 5%, 15%, 20%, as well as the 1-MGs/1,2-DGs ratio with a positive effect (from 1.7 to 2.7; Table S2, supplementary data) for 2.5, 5, 10 and 15 wt.%. The first observation is in agreement with the fact that 1,3-DGs present a more symmetrical molecular structure than 1,2-DGs, and therefore leads to better crystallization and denser packing [START_REF] Saitou | Crystallization behavior of diacylglycerol-rich oils produced from rapeseed oil. JAOCS[END_REF]. The positive effect of the 1-MGs/1,2-DGs ratio reflects the fact that MGs crystallize better and lead to more efficient gelation than DGs [START_REF] Tavernier | Relating crystallization behavior of monoacylglycerols-diacylglycerol mixtures to the strength of their crystalline network in oil[END_REF] . Monounsaturated FAs showed a negative effect for 5% and 10% and, in particular, C18:1 (-1.1, 5%; -1.3, 10%; Table S2, supplementary data), then C18:2, C16:1, and C20:1 for 10 wt.%. The effect of these unsaturated FAs is also noted by the trans FA parameter with a negative effect (-1.6, 20 wt.%;

301 Table S2, supplementary data). Whatever the chain length, the negative effect of unsaturation 302 on the crystalline packing of lipids is well known. This was clearly demonstrated, but only 303 with MGs [START_REF] Vereecken | Comparing the crystallization and polymorphic behaviour of saturated and unsaturated monoglycerides[END_REF]. Odd FAs are identified as influencing gelling 304 mechanisms. Indeed, C17:0 shows a positive effect at 10 and 20 wt.% (1.1 and 1.6

305 respectively) and C15:0, a negative effect at 15 wt.% (-1.1) (Table S2, supplementary data). A 306 similar effect is described with n-alkanes, where odd chain lengths are less organized than 307 even ones. This phenomenon, however, seems less pronounced for FAs with very long chains 308 [START_REF] Yang | Dynamic Odd-Even Effect in Liquid n-Alkanes near Their Melting Points[END_REF]. Once again, the observed effect seems to indicate a tipping point for a

309

C16 length in our mixtures, with a positive effect for C17 (> C16) and a negative one for C15

± 236792

All FAs with chain lengths less than C16:0 display significant weights with a negative influence on G'. This is valid for C8:0, C10:0, C12:0 and C14:0, and intermediate FAs,

whatever the concentrations studied. This observation confirms the tipping point previously determined at a C16 chain length.

Multiple linear regression models

To quantify the combined effects of the different structural parameters, we modeled the gelation of each AG mixture by crossing the key analytical parameters (Table S2 and Fig. S1) identified by ANNs, in MLRs, using the smallest number of non-auto-correlated key parameters. The linear equation determined by the MLR has the form as the form shown in Eq. (1).

Eq. (1). ′ % = + ∑ ×

where ′ % is the elastic modulus measured for the oil mixture at % of gelator (Table 2), are the parameters determined by least squares and are the independent analytical descriptors chosen from Table S1. Table S3 gives the values of the regression parameters, as well as the model descriptors, for the five mass concentrations investigated in this study.

From 2.5 to 15 wt.%, MLRs show that structural parameters are strictly related to C16:0, regardless of the AG class.

The C16:0 FAs have very high parameter weights for all concentrations, except for 20 wt.% (0.550, 15 wt.%; 0.921, 10 wt.%; 0.157, 5 wt.% and 0.207, 2.5 wt.%; Table S3) and a positive influence in linear regressions. C16:0-rich oleogelators (Lamemul, MDG5611, Multec Mono 9202) show a strong influence on the structuration of rapeseed oil at the lowest concentration (2.5%), with measured G' values greater than 1000 Pa (Table 2). For the highest concentration of AGs mixtures, the models do not provide a satisfactory regression with C16:0 and in this case, the C18:0 FAs show a differentiation among the AGs mixtures in linear regressions with a high weight (0.739) (Table S3). MDG0051, which has the highest rate of C18:0 (89.8%, Table S1), has the greatest influence on G' with the highest value (>7.10 6 Pa, Table 2) obtained for the highest concentration (20%). To complete these models, it is necessary to combine information on the AG classes. Thus, the 1,3-DGs/1,2-DGs ratio has a very important positive influence for concentrations of AGs mixtures of 2.5, 15 and 20 wt.% (1.080, 0.742 and 0.663, respectively) (Table S3). For other concentrations, the MGs/DGs ratio (0.548, 10 wt.%) shows an effect (Table S3) as does 1-MGs/1,2-DGs ratio (0.935, 5 wt.%) (Table S3).

MLRs show the negative influence of unsaturated FAs C18:2 (-0.514, 10 wt.%; Table S3) and C18:1 (-0.435, 5 wt.%; Table S3). MDG5611, which has a higher rate of C16:0 (59.9 g/100g)

than Lamemul (40.9 g/100g) (Table S1), but a high rate of C18:1 (20.0 g/100g) and lower G' values for all concentrations (1173 Pa, 2.5%; 3679 Pa, 5%; 155231 Pa, 10%; 503483 Pa, 15%; and 894604 Pa, 20%; Table 2). A difference can be highlighted between mono-and poly-unsaturated and the ratios of these classes. A positive effect is determined for unsaturated/saturated FA ratios (0.651, 20 wt.%) (Table S3), while a negative effect is attributed to monounsaturated/polyunsaturated FA ratios (-0.668, 20 wt.%) (Table S3). We observed that for high concentrations, unsaturated compounds do not affect gelation, because the crystalline organization is strong enough and not disturbed by the presence of these unsaturated chains. In addition, it has been observed that polyunsaturated chains disturb the organization of the AGs more than the monounsaturated ones [START_REF] Vereecken | Comparing the crystallization and polymorphic behaviour of saturated and unsaturated monoglycerides[END_REF]. The negative influence of C18:0 FAs is confirmed for low acylglyceride mixture concentrations, with low weights (-0.101 and -0.019 for 2.5 and 5 wt.%, respectively) (Table S3). The negative contribution of short-chain FAs is confirmed for intermediate concentrations (C15:0, -0.311 at 10 wt.%; C12:0 -0.209 at 15 wt.%) (Table S3). The theoretical G' values were calculated from all MLRs (Eq. ( 1) and Table S3 for the parameters and descriptors) and compared to normalized experimental values (Fig. 1). For all models, R² is greater than 0.7 for n=13, confirming regressions. The small number of AGs mixtures samples prohibits the study of interactions between descriptors or external validation, but classification of AGs mixtures remains possible on the basis of the weights studied by ANNs and confirmed by MLRs, as shown in Fig. 1. All statistical results are presented in Table S4. For the fitting validation, all MLRs give R-squared models in the 0.733-0.963 range for n=13 samples and p=4 or 5 independent analytical parameters. The p-value is below the 5% confidence interval for all MLRs, in the 2.9.10 -2 -8.6.10 -6 range. For cross-validation, R-squared for the calculated G' vs experimental regressions are in the 0.702-0.963 range and the RMS errors are in the 0.003-0.029 range. All these data provide validation of all the linear regression parameters according to the statistical criteria for our experimental space. It should be noted that the low number of commercial AGs makes it impossible to evaluate external cross-validation using different test and regression datasets. These models can therefore only be used for the analysis of the interactions of the variables corresponding to the analytical parameters [START_REF] Gramatica | On the development and validation of QSAR models[END_REF][START_REF] Konovalov | Robust cross-validation of linear regression QSAR models[END_REF].

Principal component analysis

The principal component analysis (PCA) of rapeseed oil samples gelled with low concentrations of each AGs mixture (2.5 wt.%; Fig. 2) shows a group of 4 oleogelators (Lamemul: AG3; glycerol monostearate: AG7; MDG 5611: AG8 and Multec Mono 9202: AG13) in the same space, defined by the 3 components PC1, PC2 and PC3 representing 90% of the observed variance, using C16:0, the ratio 1,3-DGs/1,2-DGs and the ratio 1,2-DGs/TGs as the main descriptors favorable to interactions structuring the oil matrix. The positive effect of the 1,2-DGs/TGs ratio reflects the better gelling properties of 1,2-DGs compared to triglycerides, even if its crystalline organization is less dense than that of 1,3-DGs [START_REF] Tavernier | Relating crystallization behavior of monoacylglycerols-diacylglycerol mixtures to the strength of their crystalline network in oil[END_REF]. This group shows that a mixture of AGs containing more than 32 wt.% of C16:0, a 1,3-DGs/1,2-DGs ratio greater than 32.510, and a 1,2-DGs/TGs ratio greater than 2.510 explain the gelling properties at low concentrations. Longer chain FAs (i.e., C18:0) do not explain classification, as the negative contribution of this descriptor indicated in Table S3 (a4, 2.5%). Thus, MDG 0091 (AG9), MDG 0096 (AG4) and MDG 0051 (AG11), with the highest concentrations in C18:0 FAs (ranging from 78.5% to 90.61%) and lower values for the C16:0 parameters and the 1,3-DGs/1,2-DGs ratio, are classified in the group of weak oleogelators (Fig. 2, group B). If the MGs correspond to the major species for most of the AGs studied (Table S1, c), Supplementary data), it is clear that DGs can contribute to the classification.

The 1,3-DGs/1,2-DGs ratio shows that the balance between these two isomers is an important factor to be considered for effective oleogelating properties at low concentrations. The PCA analysis for the highest concentration of alkyl glycerides (20 wt.%) (Fig. 3.)

shows a group of 5 oleogelators (Lamemul: AG3; Glycerol monostearate: AG7; MDG 0298: AG10; MDG 0051: AG11 and Multec Mono 9202: AG13) in the same space, defined by the 3 components PC1, PC3 and PC4 representing 70% of the variance observed. These 3 PCs are the product of the key parameters identified with MLRs. Among these parameters, C18:0, the unsaturated/saturated FA ratio and the 1,3-DGs/1,2-DGs ratio show a positive weight in the interactions favorable to G' and obtaining strong oleogels (Table S3). But the negative contribution of the monounsaturated/polyunsaturated FA ratio and the 2-MGs/TGs ratio must be taken into consideration to understand the interactions that explain the gelation properties.

In fact, MDG 0051 (AG11, Fig. 3), which gives the highest value of G' for a 20 wt.% concentration in rapeseed oil (7.74.10 6 Pa), shows a very high rate of C18:0 (89.9%) as well as a high 1,3-DGs/1,2-DGs ratio (3.38), but MDG 0091 (AG9, Fig. 3) with the highest rate of C18:0 (90.6%) and a high 1,3-DGs/1,2-DGs ratio (2.94) is classified in Group C with an intermediate G' value (2.01.10 5 Pa). The negative contribution of the 2-MGs/TGs ratio for MDG 0091 (67.0) is much higher than that of MDG 0051 (0.87). It is important to note that the contribution of C16:0 no longer explains the properties of each AG for a 20 wt.% concentration. It appears from this classification that C18:0 is not solely responsible for the good properties of the oleogelators. We also observed that the high weight given to the unsaturated/saturated FA ratio explains the behavior of MDG 0298 (AG10, Fig. 3). The latter shows a very low level of C18:0 (2.7%), but the highest unsaturated/saturated FA ratio (11), related to a very high unsaturated C18:1 FA level (83.0%). Moreover, this oleogelator shows a very low 2-MGs/TGs ratio (3.81), which gives less of a negative contribution to interactions. Therefore, it appears that AGs containing C18:1 long-chain FAs, mainly in the form of MGs, may be considered as efficient oleogelators valuable at high concentrations, but gelling mechanisms with these same species in the mixture at low concentration (2.5 wt.%) are not possible (Fig. 2 (B),AG10). This identified high concentration complexity is reflected in the appearance of the other two G' intermediate spaces (C and D, Fig. 3). Space C contains 3 oleogelators MDG 0096: AG4; MDG0091: AG9 and MDG 5616: AG12, which have high levels of C18:0, but the three of them have the highest monounsaturated/polyunsaturated FAs or 2-MGs/TGs ratios of all the AGs studied. The PCA study of the other 3 concentrations: 5, 10 and 15 wt.% (Fig. 4) shows continuity with the study of the 2.5 wt.% concentration. The classification is obtained in two groups according to the gelling properties and the G' values measured (Fig. 4). These 3 concentrations are classified according to the components PC1, PC2, PC3 representing 90% of the variance for 5 wt.%; PC1, PC2 and PC5 representing 74% of the variance for 10 wt.% and PC1, PC2 and PC4 representing 86% of the variance for 15 wt.%. For each of the MLRs and principal components, C16:0 appears to have a strong positive influence on the interactions contributing to G', associated with the 1,3-DGs/1,2-DGs ratio for 5 and 15 wt.% and the MGs/DGs ratio for 10 wt.%. For these intermediate concentrations, the influence of C18:1 and C18:2 FAs appears negatively in interactions contributing to G'.

Conclusion

This study shows that the classification of commercial MGs mixtures is possible using We have demonstrated that by combining 4 or 5 analytical parameters it is possible to classify different commercial MGs mixtures according to their gelling properties in vegetable oil and according to their concentration. This classification showed that 3 products were grouped and can be considered as strong oleogelators for all the concentrations tested: Lamemul: AG3; glycerol monostearate: AG7 and Multec Mono 9202: AG13.

This study demonstrates that the combination of a neural network study associated with multiple linear regression models and a principal component analysis constitutes a powerful tool for the study of complex mixtures. Such an approach could be applied to a large number of fields (food, cosmetics or pharmaceuticals) for industrial applications but also for a more fundamental comprehension of multiple interactions in complex mixtures.

  to 10 C); long-chain FAs (12 to 22 C) and very long-chain FAs (> 22 C). 1 H NMR data were introduced for the determination of the amount of the different AGs and their isomers: TGs; 1,2-DGs; 1,3-DGs; 1-MGs; 2-MGs. Specific ratios were calculated from the raw GC data: unsaturated/saturated FAs; monounsaturated/saturated FAs; monounsaturated/polyunsaturated FAs; FAs polyunsaturated/saturated; long-chain/intermediate FAs; very long/long-chain FAs; long chain/intermediate FAs; intermediate/long-chain FAs; intermediate/very long-chain FAs;

Fig. 1 .

 1 Fig. 1. Calculated vs experimental G' for the 13 AGs mixtures by varying mass concentrations from 20 to 2.5 wt.% in rapeseed oil (from A to E, respectively).

Fig. 2 .

 2 Fig. 2. PCA based on the key analytical parameters for 2.5 wt.% of AGs in rapeseed oil clustered in high measured G' values > 8.10 2 Pa (A) and low G' < 8.10 2 Pa (B).

Fig. 3 .

 3 Fig. 3. PCA based on the key analytical parameters for 20 wt.% of AGs in rapeseed oil clustered in high measured G' values > 3.10 6 Pa (A), low G' < 4.10 3 Pa (B), and intermediate G' from 2.10 5 to 5.10 5 Pa (C) and from 8.5.10 5 to 9.5.10 5 Pa (D).

Fig. 4 .

 4 Fig. 4. PCA based on the key analytical parameters for 5 wt.% of AGs in rapeseed oil clustered in high measured G' values > 3.610 3 Pa (A), low G' < 2.4.10 3 Pa (B), 10 wt.% clustered in high measured G' values > 1.5.10 5 Pa (A), low G' < 7.4.10 3 Pa (B), and 15 wt.% clustered in high measured G' values > 5.10 5 Pa (A), low G' < 2.4.10 5 Pa (B).

  analytical descriptors based on 1 H NMR and GC. The use of ANNs made it possible to identify the key analytical parameters that influence gelling properties depending mainly on the nature of the oleogelator involved and their respective concentrations (from 2.5 to 20 wt.%). The oleogelators studied are mainly composed of MGs, among them 1-MGs.The analytical descriptors used in this study are adapted to this class of lipids and allow classification without studying the possible synergies induced by the variations in the minor species, such as the different positions of FAs on DGs or TGs. MLR models have confirmed the weight and influence of the selected parameters on the gelation and hardness of the final dispersion. The QSPR study showed that MGs-rich AGs with C16:0 FAs, high 1,3-DGs/1,2-DGs ratios and low intermediate AGs (<C15 FAs) present good gelling properties at low concentrations up to 15 wt.% in rapeseed oil. The MGs or DGs containing C18:0 present in these mixtures were not able to contribute positively to gelation (G') at low concentrations, but were very effective at high concentrations, especially at 20 wt.%. It appears that these same species contribute negatively to G' at low concentrations. This study also showed that MGs containing C18:1 strongly contribute to efficient gelation of the oil at high concentrations, but not at the lowest concentrations, where they contribute negatively to G'.

Table 1 .

 1 Macroscopic observations of the physical state of oily mixtures at room 243 temperature after 24 hours, made with a growing amount of commercial AGs mixtures. No unsaturated FAs are found and chain lengths are exclusively limited to C8:0 and C10:0. The other commercial AGs mixtures behave as gelators, showing variable gelation properties depending on their composition and concentration. We observed two categories of products: Lamemul, MDG 0096 and Multec Mono 9202, which show high gelation properties starting from 2.5 wt.%, and the other eight products which start to gelate

	244
					Concentrations of AGs mixtures (wt. %)
		AGs mixtures	Sample
				2.5	5	10	15	20
			code
	C8-10 MGs intermediate Type II AG1	0	0	0	0
		C8-10 MGs intermediate Type I AG2	0	0	0	0
		Lamemul	AG3	1	1	1	1
		MDG 0096	AG4	1	1	1	1
		Glyceryl dibehenate	AG5	0	0.5	1	1
		HS-glycerides	AG6	0	0	0.5	0.5
		Glycerol monostearate	AG7	0	0.5	1	1
		MDG 5611	AG8	0.5	1	1	1
		MDG 0091	AG9	0	1	1	1
		MDG 0298	AG10	0	0	0.5	1
		MDG 0051	AG11	0	0.5	1	1
		MDG 5616	AG12	0	0.5	1	1
		Multec Mono 9202	AG13	1	1	1	1
	245	0: liquid state; 1: gel state; 0.5: intermediate state, liquid+gel
	246
	247	C8-10 MGs intermediate Type I and Type II were selected for their absence of
		gelation, confirmed by the liquid state observed for all concentration in Table 1. The main

248

MGs of C8-10 Type I and Type II are 1-MGs (64.4% and 84.4%, respectively, Table

S1

, 249 supplementary data). which allow a liquid oil phase to appear if the tubes are inverted. This is the case for Glyceryl dibehenate, HS-glycerides, Glycerol monostearate, MDG 5611, MDG 0298, MDG 0051 and MDG 5616 for concentrations below 20% (Table

Acknowledgments

The authors gratefully acknowledge the Nouvelle-Aquitaine Regional Council, the Délégation Régionale à la Recherche et à la Technologie (DRRT Nouvelle Aquitaine) and the Agence Nationale pour la Recherche Technologique (ANRT) for their financial support. The authors would also like to thank Franck Dejean from ITERG (Pessac, France) for his helpful support on the GC analysis.


Joris Fournel
email: jorisfournell@gmail.com

Axel Bartoli

David Bendahan

Maxime Guye

Monique Bernard

Elisa Rauseo

Mohammed Y Khanji

Steffen E Petersen

Alexis Jacquier

Badih Ghattas

Medical image segmentation automatic quality control: a multi-dimensional approach

Keywords: Medical Image Segmentation Automatic Quality Control, Multi-dimensional quality control, CMR image segmentation, Deep Learning

In clinical applications, using erroneous segmentations of medical images can have dramatic consequences. Current approaches dedicated to medical image segmentation automatic quality control do not predict segmentation quality at slice-level (2D), resulting in sub-optimal evaluations. Our 2D-based deep learning method simultaneously performs quality control at 2D-level and 3D-level for cardiovascular MR image segmentations. We compared it with 3D approaches by training both on 36540 (2D) / 3842 (3D) samples to predict Dice Similarity Coefficients (DSC) for 4 different structures from the left ventricle, i.e., trabeculations (LVT), myocardium (LVM), papillary muscles (LVPM) and blood (LVC). addReviewerThe 2D-based method outperformed the 3D method. At the 2Dlevel, the mean absolute errors (MAEs) of the DSC predictions for 3823 samples, were 0.02, 0.02, 0.05 and 0.02 for LVM, LVC, LVT and LVPM, respectively. At the 3D-level, for 402 samples, the corresponding MAEs were 0.02, 0.01, 0.02 and 0.04. The method was validated in a clinical practice evaluation against semiqualitative scores provided by expert cardiologists for 1016 subjects of the UK BioBank. Finally, we provided evidence that a multi-level QC could be used to enhance clinical measurements derived from image segmentations.

Introduction

Manual segmentation of medical images is time consuming and subject to inherent between-and within-operator variabilities. Artificial intelligence (AI)-based methods control (QC) of segmentations is thus critical in order to validate the robustness of the corresponding methods.

On that basis, Medical Image Segmentation Automatic Quality Control (MISAQC) methods have been growingly developed over the last few years. MISAQC methods intend to automatically evaluate the quality of a given segmentation and to track and potentially discard the segmentations which do not fulfill the quality criteria.

The quality of AI-based segmentation can be assessed using various metrics such as Dice Similarity Coefficient (DSC), Hausdorff distance, volume/mass error, etc ... This quality is considered to be good if these indices lie within a commonly acceptable range in the corresponding field. The computation of those metrics is based on a comparative analysis between the predicted segmentation and the ground truth manual segmentation. On the contrary, MIS-AQC methods aim at predicting a quality metric (e.g. DSC) for a given image and its predicted segmentation, and thus do not use the ground truth manual segmentation.

So far, very few approaches have been reported regarding MISAQC ( [START_REF] Albà | Automatic initialization and quality control of large-scale cardiac mri segmentations[END_REF][START_REF] Audelan | Unsupervised quality control of image segmentation based on bayesian learning[END_REF][START_REF] Kohlberger | Evaluating segmentation error with-out ground truth[END_REF][START_REF] Valindria | Reverse classification accuracy: Predicting segmentation performance in the absence of ground truth[END_REF]). A regression model was used to assess the quality for automated segmentation of lung and liver images [START_REF] Kohlberger | Evaluating segmentation error with-out ground truth[END_REF]. The DSC was used as a QC index and variable correlations between the predicted and real DSCs for the lung (r=0.85) and the liver (r=0.54) have been reported. [START_REF] Audelan | Unsupervised quality control of image segmentation based on bayesian learning[END_REF]) used an unsupervised bayesian approach and reported good correlations between the real and predicted DSC scores for brain tumors (r=0.69) and the left ventricle myocardium (r=0.78).

Furthermore, [START_REF] Robinson | Automated quality control in image segmentation: Application to the uk biobank cardiovascular magnetic resonance imaging study[END_REF] showed that the reverse classification accuracy (RCA, [START_REF] Valindria | Reverse classification accuracy: Predicting segmentation performance in the absence of ground truth[END_REF]) approach could provide good DSC predictions without the need of a large annotated dataset. However, the time required for the RCA analysis of a single segmentation (11 min) was considered prohibitive for real-time applications [START_REF] Robinson | Real-time prediction of segmentation quality[END_REF]. For the segmentation of cardiac magnetic resonance images, [START_REF] Robinson | Real-time prediction of segmentation quality[END_REF] used a 3D convolutional neural network (CNN) in order to predict the DSC values of 3D segmentations. The corresponding results were very promising in terms of predictive accuracy, with mean absolute errors (MAEs) between the real and predicted DSCs of 0.03 ± 0.04 for the whole heart and a very short processing time (less than a second).

However, the segmentation quality assessment was based on the 3D DSC and this approach can actually be questioned. The 3D DSC is a volume-related information that is very useful to evaluate the global quality of a segmentation but it does not contain any specific information related to the individual segmented slices which are part of the volume of interest. In other words, errors cannot be accurately localized using a 3D DSC-based method. Significant segmentation errors can occur at the 2D level and be fully ignored at a 3D level thereby overestimating the segmentation quality. The localization of segmenta- tion errors at a slice-level could offer the opportunity for corrective processes thereby providing more accurate volumetric measurements and potentially improved diagnosis. On that basis, we hypothesized that a MISAQC tool would largely benefit from both 3D and 2D evaluations.

In the present paper, our goal was to present a new MISAQC tool having the following properties:

• Easily trainable using only 2D CNNs.

• Provides quality assessment at 2D-level, thus localizes slices where the segmentation is erroneous.

• Helps to automatically correct clinical measurements biases induced by erroneous segmentations.

• Provides quality assessment at 3D-level from a mathematical combination of 2D-level predictions.

Materials and Methods

Study design

This multicenter retrospective study was approved by the local institutional review board (N °IRB CRM -1907-02è) in accordance with the guidelines outlined in the Declaration of Helsinki.

A CNN-based method was developed for the automated assessment of segmentation quality for the LV in the enddiastolic (ED) phase from CMR images [START_REF] Bartoli | Deep learning-based automated segmentation of the left ventricular trabeculations and myocardium on cardiac mr images: A feasibility study[END_REF].

The initial dataset included 4290 images from 449 different subjects. Patients were classified as healthy, dilated, hypertrophic, hypertrabeculation by the most experienced investigator (A.J. -20 years of experience). Healthy refers to a normal cardiovascular examinations. Dilated indicated a dilated cardiomyopathy (DCM) as previoulsy described [START_REF] Japp | The diagnosis and evaluation of dilated cardiomyopathy[END_REF]. Hypertrophic refers to hypertrophic cardiomyopathy (HCM) which was defined based on the myocardial thickness on the ED short-axis [START_REF] Elliott | 2014 esc guidelines on diagnosis and management of hypertrophic cardiomyopathy: The task force for the diagnosis and management of hypertrophic cardiomyopathy of the european society of cardiology (esc)[END_REF]. Excessive trabeculation cardiomyopathy (ETCM) was defined according to [START_REF] Petersen | Left ventricular non-compaction: insights from cardiovascular magnetic resonance imaging[END_REF] criteria i.e. a double-layer myocardium aspect and a noncompacted layer to compacted layer thickness > 2.3 on ED long-axis. In our dataset, 259 subjects were healthy subjects, 39 were DCM, 92 were HCM and 59 were ETCM. In Table 1 are displayed the patient characteristics for this initial dataset.

CMR Data

The CMR examinations were performed at 1.5T using 2 different scanners i.e. Ingenia ® (Philips Health System, Best, the Netherlands) and Avanto ® (Siemens Healthcare, Erlangen, Germany). Patients were positioned supine in the scanner with a multi-channel body array coil positioned on the top and a spine array coil positioned on the bottom. MR acquisition was gated to ECG and occurred during an inspiratory breath hold. A balanced turbo fieldecho sequence and a balanced steady-state-free precession sequence was used with the Ingenia and Avanto scanners respectively. In both cases, images were acquired in the short-axis view in order to cover the LV from base to apex. The following parameters were used : in-plane resolution : 1.5 mm 2 (Ingenia) and 2.35 mm 2 ,(Avanto) slice thickness : 7-8 mm, gap between slices :10 mm, flip angle : 30. The short-axis image stack consisted in 8 to 17 slices depending on the scanner, patient height, cardiac anatomy and morphology. The ED frame at each imaging level was retained for further image analysis. Data of each patient was de-identified before the analysis.

Manual segmentation and reference measures

Manual image segmentation was undertaken by a trained observer (A.B., who has 5 years of experience) using a post-processing software previously validated by Bricq et al. (Bentatou et al., 2018;[START_REF] Bricq | Semiautomatic detection of myocardial contours in order to investigate normal values of the left ventricular trabeculated mass using mri[END_REF][START_REF] Frandon | Semi-automatic detection of myocardial trabeculation using cardiovascular magnetic resonance: Correlation with histology and reproducibility in a mouse model of non-compaction[END_REF]. The LV structures were manually segmented so as to obtain four labels i.e: blood cavity (LVC), myocardium (LVM), papillary muscles (LVPM) and trabeculations (LVT). The non segmented part of each image was assigned to a fifth label: background (BG). The corresponding results (Manual) were used as the reference standard.

Methods

The aim of a MISAQC approach is to learn to predict the quality of a segmentation. In order to train a supervised model for that purpose, we first need to build a dataset containing couples of images/segmentations with varying qualities for the segmentations as shown in figure 2. Once this dataset is built, different CNN models are trained and compared for the 2D and 3D DSC prediction. Our approach (illustrated in figure 3) consists of two steps. First, we predict 2D indices of quality for each segmented slice. Second, we use the 3D DSC formula to estimate its components using the 2D indices predicted in the first step. We compare our approach to a direct 3D approach.

MISAQC dataset construction

To generate the MISAQC dataset, we used an approach similar to the one proposed by [START_REF] Robinson | Real-time prediction of segmentation quality[END_REF] illustrated in figure 2. Twelve versions of a same CNN model (see [START_REF] Bartoli | Deep learning-based automated segmentation of the left ventricular trabeculations and myocardium on cardiac mr images: A feasibility study[END_REF] for details about the architecture) differing by the number of feature maps per convolutional layer, the proportion of the original dataset used for training, and the number of training epochs are implemented. The three hyperparameters of the CNNs were selected randomly, resulting in models producing segmentations of variable qualities. To train the CNNs we used the 4290 2D CMR images for which manually segmented masks were available. Each of these models was then used in inference to generate segmentations for the whole original dataset. The corresponding 2D and 3D DSC scores were computed for all the segmentations. Thus for each original 2D image, twelve segmentations are available together with their corresponding DSC. This resulted in a total of 5388 3D segmentations, corresponding to 51480 2D segmentations.

The next step guarantees that the obtained dataset is balanced with respect to the quality of the segmentations. Therefore, we sample from this artificial dataset triplets (image, mask, DSC) stratifying over LVT 3D DSC distribution, using the following bins: [0, 0.2], [0.2, 0.3], [0.3, For each segmentation, 5 one-hot-encoded masks were generated, one for each class (BG, LVM, LVC, LVT and LVMP), and concatenated with the CMR input. Segmentations and CMR inputs were previously cropped around the LV (see [START_REF] Bartoli | Deep learning-based automated segmentation of the left ventricular trabeculations and myocardium on cardiac mr images: A feasibility study[END_REF] for more details) and the CMR inputs were pre-processed with Contrast Limited Adaptive Histogram Equalization (CLAHE). The 2D inputs size was 128x128x6 while the 3D inputs were resized to a 128x128x8x6 using nearest neighbour interpolation in order to account for the variable number of slices per subject. The 2D inputs are used for our method and the 3D inputs for the other approaches for comparative purposes. This process resulted in a total of 36540 2D and 3840 3D samples for training and validation. The testing phase was performed on 3823 2D and 402 3D samples.

3D DSC as a function of 2D indices

We will show how the 3D DSC may be recovered from the 2D DSC and the mean volume similarity fraction (MVSF). Let Y and Z be two 3D regions (Y the groundtruth 3D segmentation and Z the corresponding automated segmentation) composed of n 2D stacks, or equivalently:

Y = (y i ) i∈ 1,n
and:

Z = (z i ) i∈ 1,n
where y i and z i are 2D regions. The 3D DSC and MVSF between Y and Z are defined as: Similarly, let DSC(y i ,z i ) and MVSF(y i ,z i ) be the DSC and MVSF between the 2D slices y i and z i . DSC(Y ,Z) is used as the quality index of the subject-level automated segmentation Z, while DSC(y i ,z i ) serves the same purpose for the slice-level segmentation z i . The goal is to predict these metrics in the absence of the ground truth, Y .

DSC(Y, Z) = 2|Y ∩ Z| |Y | + |Z| (1) M V SF (Y, Z) = 2(|Y | -|Z|) |Y | + |Z| (2)
Each example from the learning sample was composed of an MR image (x i ), the corresponding automated (z i ) and ground-truth segmentation (y i ). In order to evaluate the quality of Z, we wish to predict both scores (DSC and MSVF) relative to each pair (y i , z i ) i∈ 1,n .

A model, denoted by M , is trained such that:

M (x i , z i ) = [ DSC(y i , z i ), M V SF (y i , z i )]
This model can be used to predict, for each new pair consisting of an MR image and its automated segmentation, a value for both scores without requiring the ground-truth segmentation. We will now see how the recovery of the 3D DSC can be achieved. Considering the expression for the 2D DSC:

DSC(y i , z i ) = 2|y i ∩ z i | |y i | + |z i | (3)
we obtain the following estimation:

|y i ∩ z i | = 1 2 DSC(y i , z i )(( |y i | + |z i |) (4)
And considering the 2D MVSF expression:

M V SF (y i , z i ) = 2(|y i | -|z i |) |y i | + |z i | (5)
we can estimate |y i |:

|y i | = |z i | 2 + M V SF (y i , z i ) 2 -M V SF (y i , z i ) (6)
The 3D DSC predictions can be recovered as follows:

DSC(Y, Z) = 2 |Y ∩ Z| |Y | + |Z| (7)
Where

|Y | = i∈ 1,n |y i |, |Z| = i∈ 1,n |z i | and |Y ∩ Z| = i∈ 1,n |y i ∩ z i |.

CNN architectures for DSC prediction

To implement our approach, we suggest a CNN architecture for the prediction of both 2D DSC and MVSF. We use then the formulas from the previous section to obtain the 3D DSC values. In order to show that the performance of our approach is not dependant on a specific architecture, two types of backbone architectures were used i.e QC ResNet and QC U-Net (figure 4). The rationale behind the QC U-Net backbone type was to build an architecture that would implicitly reproduce the steps performed by a human expert to produce a DSC for a given segmentation and related medical image. The expert would (1) segment the image manually and then (2) compute the DSC between his segmentation (considered as ground truth) and the segmentation to evaluate. Hence, we considered the image segmentation task as a necessary step in the process of a DSC computation task for a human expert. Having this in mind, the QC U-Net backbone architecture was designed as a two-stages process, i.e., a U-Net encoder-decoder block (originally designed for image segmentation) augmented with two layers aiming at directly predicting the DSC values. Based on empirical preliminary measurements, the QC ResNet version was designed with a fewer number of features maps than in the original ResNet-50 [START_REF] He | Deep residual learning for image recognition[END_REF]. Finally, to compare our approach to methods directly predicting the 3D DSC [START_REF] Robinson | Real-time prediction of segmentation quality[END_REF], we used similar networks adapted for 3D data structure. The 2D-based method was identified as 2D R , where R stands for reconstruction given that the 3D DSCs were reconstructed from 2D indices. For each class, the 2D R method performed a multi-level quality control: it predicted a 2D DSC for each segmented slice together with a reconstructed 3D DSC for the whole volume. For each implementation of our 2D R method, 2D images and the corresponding segmentations were taken as inputs while a 2D CNN was trained to predict both the 2D DSC and the 2D MVSF scores for each class (BG, LVM, LVC, LVT and LVPM). On that basis, the last layer was composed of ten units: five for the 2D DSC scores and five for the 2D MVSF scores. For the 3D direct approach, taking 3D images and their segmentations as inputs, networks were trained to predict the five 3D DSC resulting in five units for the last layer. As indicated in figure 4, the remaining differences between the 2D and 3D networks for a given backbone architecture type were the kernel size of the convolutional layers, the pooling size in pooling layers and the input shape.

Implementation details

Each 2D network was trained for 80 epochs (until no significant gain was observed on the validation set) using Adam optimizer at a learning rate of 0.001. Similarly, each 3D network was trained for 50 epochs using Adam optimizer at a learning rate of 0.001. For each network the loss function was the mean squared error (MSE). No data augmentation was used.

Evaluation

Slice-level evaluation

The performances of the 2D R method implemented with each backbone (QC ResNet and QC U-Net) were initially assessed at the slice-level using 3823 segmented 2D images of the test dataset. For each class, the evaluation was carried out on the basis of the mean absolute error (MAE) which was calculated between the predicted 2D DSCs and the ground-truth values. For each class, the accuracy indicated the discriminative capacity of a given method to distinguish "bad quality" from "good quality" 2D segmentations based on the predicted 2D DSC values. To compute the accuracy, predicted and observed DSC values were encoded as "bad quality" or "good quality" according to the following thresholds: 0.7 for LVM and LVC according to [START_REF] Robinson | Real-time prediction of segmentation quality[END_REF]. Considering that LVPM and LVT (as illustrated in figure 5) are prone to lower DSC values [START_REF] Bartoli | Deep learning-based automated segmentation of the left ventricular trabeculations and myocardium on cardiac mr images: A feasibility study[END_REF] (as illustrated in figure 5), a lower threshold was used, i.e., 0.45 for LVPM and 0.35 for LVT. Pearson correlation coefficients were also computed in order to further evaluate the performances at the slice-level.

Subject-level evaluation

The performance of the 2D R method was similarly evaluated for each backbone at a 3D-level using the 402 segmented 3D images of the testing dataset. The correspond- ing results were compared to those obtained with a 3D direct approach.

Clinical practice evaluation

Our trained 2D R QC U-Net model was communicated to the research team of Steffen Petersen that used it to control the quality of 1016 segmented subjects from the UK BioBank. For this dataset, trained cardiologists also assigned a semi-qualitative quality control score (QC 1=Good, QC 2=Sub-optimal but still usable, QC 3=Poor) for each class and slice. At the 3D level, QC 2 and QC 3 are merged and patient is classified as QC 2/3 if at least one slice belongs to either QC 2 or QC 3. The agreement between these human scores and our model predictions was assessed.

Statistical tests

Statistical comparison were performed (for mean differences) using Wilcoxon signed-rank tests for paired samples, unpaired t-test for unpaired samples. Two-sample Kolmogorov-Smirnov test was used for testing whether two independent samples are drawn from the same continuous distribution. Differences were considered as significant for p-values lower than 0.05.

Results

Slice-level results

The comparative analysis between the different backbones is summarized in Table 2. As illustrated in Table 2, regardless of the backbone and for the whole set of classes but LVT, the MAEs computed at the slice-level i.e. 2D were systematically lower than 0.023. The MAE related to the LVT class was slightly larger i.e. 0.047.

For the QC U-Net backbone, the MAE associated with the 2D R method was 0.005 for the BG 2D DSC prediction, 0.023 for LVM, 0.019 for LVC, 0.047 for LVT and 0.020 for LVPM. Of note the largest error was found for a small and scattered class i.e. LVT.

When comparing backbones, MAEs values were found significantly smaller for the QC U-Net 2D R for all the classes but LVT.

Accuracies in distinguishing bad from good quality slicelevel segmentations were larger than 96.1% for 6 out of the 8 tested classes and larger than 92.7% for the LVT class. Overall, the segmentation quality predictions at the slice-level were characterized by low MAEs values regardless of the backbone and the class.

Subject-level results

As illustrated in Table 3 andfigure 6, the MAEs computed from the ground-truth and the predicted DSC val-ues were systematically and significantly lower for the 2D R variant and so regardless of the backbone used. A high level of predictive accuracy was consistently achieved for the whole set of classes. Regarding the QC ResNet backbone, the 2D R (3D) method provided a MAE of 0.004 (0.011) for the BG 3D DSC prediction, 0.016 (0.028) for LVM, 0.011 (0.016) for LVC, 0.022 (0.022) for LVT and Regardless of the backbone, errors of the 2D R method were at least twice smaller than those from the 3D direct approach and so for the whole set of classes but LVPM. The results of the corresponding statistical analysis are summarized in (Table 3). For the QC ResNet backbone, all the MAEs values were significantly smaller for the 2D R version. For the QC U-Net backbone it was the case for the whole set of class but LVPM. Accuracy values and pearson correlation coefficients for each backbone and method are summarized in Table 4. Considering the QC ResNet backbone, the accuracy values ranged from 90.7% (LVT) to 96.5% (LVPM) for the 3D method, while they ranged from 97.0% (LVC) to 97.7% (LVM) for the 2D R method. The correlation coefficients ranged from 0.955 (LVT) to 0.983 (LVM) for the 3D method while they ranged from 0.965 (LVPM) to 0.996 (LVC) for the 2D R method. These values were systematically superior for the 2D R version as compared to the 3D (Table 4). As a matter of example, for the LVM class the accuracy related to the 2D R method was 97.7% i.e. 5.0% larger than the corresponding 3D method's value (92.7%).

Similar results were found for the QC U-Net backbone and for the whole set of structures but LVPM (Table 4). The accuracy values from the two implementations of the 2D R method were larger than 96.5% regardless of the evaluated class while 6 out of 8 values were lower than 96.5% for the two implementations of the 3D direct approach.

MAEs at the 3D-level were generally inferior to the corresponding values at the slice-level for the 2D R method. As an example, the MAE value for the LVC class was 0.012 at the 3D-level for QC U-Net 2D R and 0.019 at the 2D-level. Of note, this was not the case for the 3D direct approach for which MAEs at the 3D-level were generally inferior than those at the 2D-level. As an illustration, the MAE for the LVM (LVC) class was 0.037 (0.029) at the 3D-level for QC U-Net 3D whereas it was 0.023 (0.024) at the 2D-level for QC U-Net 2D R .

Scatter plots, fitted regression line and R scores of the true 3D DSC against the predicted values from the QC U-Net 2D R are illustrated in Figure 7. The scatter plots illustrate the capacity of the method to predict the 3D DSC with similar accuracy in the low, middle and high range DSC values. R scores are significantly high for all classes (from 0.932 for LVPM to 0.989 for LVC).

Table 5 displays the 3D DSC MAE on the test dataset for the QC U-Net 2D R , for each patient class (DCM, HCM, ETCM, Healthy, Non-Healthy). The Non-Healthy group was composed of the patients that were not in the Healthy group. The predictions for the Healthy group generally presented smaller MAE values, but the differences were only significant for the LVC (0.010 for Healthy subjects against 0.016 for Non-Healthy, p=7.51e-4) and LVPM (0.034 against 0.049, p=8.64e-3). Overall, the method presented low MAE values for all groups and classes.

Table 6 gathers the mean processing time of each method at both 3D and 2D-levels for the 2D R version. The 2D R method had a 300 ms inference time while the 3D direct approach took about 100ms per inference.

Finally, in order to assess the impact of the training sample size, our 2D R QC-UNet approach was trained with 10%, 33%, 66% and 100% of the total training dataset. The corresponding 3D DSC MAE for LVM, LVC, LVT and LVPM are illustrated in figure 8. One can observe a substantial 3D DSC MAE reduction when the learning dataset size increases.

Clinical practice results

Figure 9 overlays the density plots of the predicted 2D DSCs for the three semi-qualitative groups (QC 1, QC 2 and QC 3) provided by the expert for each class.

First, the distribution of the predicted 2D DSCs for QC 1 always significantly (according to Kolmogorov-Smirnov test, p-values < 0.05) differed from those of the lower quality groups (QC 2 and QC 3) and so for all classes. This was especially the case for LVM and LVC which are classes that are generally present in all slices contrarily to LVT and LVPM. When a class is absent in a segmented slice, the predicted 2D DSC is zero for that slice: this explains the small mode in zero appearing for QC 1, much more present for LVT and LVPM. Secondly, the shift between the distributions of the different groups was in the same order as the scores given by the human expert (QC 3 at the left, QC 2 in the middle and QC 1 at the right in the DSC range of values) except for the LVT class for which the QC 2 and QC 3 groups had a comparable distribution shape. Moreover, all the distributions significantly differed (QC 1 vs QC 2, QC 1 vs QC 3 and QC 2 vs QC 3) except for LVT and LVPM (for QC 2 vs QC 3) where p-values were of 0.99 and 0.64 respectively.

Figure 10 shows the corresponding density plots of the predicted 3D DSC for each class. As expected, the differences between the distributions were less important than those of 2D DSCs. However, regardless of the class, the 3D DSC distribution of the QC 2/3 group was always found to be significantly different, more shifted to the left and with larger standard deviation compared to the distribution of the QC 1 group.

Automatic correction of segmentation model failure

We conducted an additional analysis in order to understand how erroneous segmentations in the midslice section could affect the 3D DSC values and the derived volumetric measurements. A simple approach is then suggested to correct the measurements in this case.

We consider the scenario of a single short-axis slice segmentation with a very poor quality in the mid-slice section. In this particular scenario, the 3D DSC value might be acceptable whereas the 3D volumetric measurement would not. In fact, measurements derived from mid-slices, where the largest section of the LV is commonly located, are likely to contribute the most to the overall volumetric quantification.

In an attempt to address this issue, we used of the CNNpredicted segmentations for the test dataset reported in [START_REF] Bartoli | Deep learning-based automated segmentation of the left ventricular trabeculations and myocardium on cardiac mr images: A feasibility study[END_REF], which was composed of 150 subjects. These segmentations were of high quality and the corresponding clinical quality metrics were accessible in [START_REF] Bartoli | Deep learning-based automated segmentation of the left ventricular trabeculations and myocardium on cardiac mr images: A feasibility study[END_REF]. For each subject, the original 3D predicted segmentation was corrupted based on a random replacement of the mid-slice segmentation by another segmentation from our MISAQC dataset. As a result, we put together two sets of 150 3D segmentations : the original dataset without deterioration (w/o D), and the dataset with deteriorations (D), for which the mid-slices were systematically erroneous. The comparative analysis between the two datasets in terms of 3D DSC values and relative volume errors (RVE) allowed to quantify the effect of midslice segmentation errors on those two metrics and to check if the metrics were divergently affected.

Furthermore, this study provided an opportunity to investigate whether our 2D-based QC model could be used to correct the overall volume measurement. In order to do so, we used our 2D R QC U-Net model to predict the 2D DSC and 3D DSC for every segmented slice and subject, respectively. Each segmented slice was classified as "good" or "wrong" based on the predicted 2D DSC: segmented slices with a predicted 2D DSC lower than 0.3 for the BG or LVM or LVC class were classified as "wrong". We chose these classes as they were most likely present in all slices as opposed to the LVT and LVPM classes. The most basal and apical slices were also excluded. We then replaced the associated volume measurements in the segmented slices identified as "wrong" by the mean of measurements in the two adjacent slices identified as "good" and so for the whole set of classes. The corresponding measurements were referred to as deteriorated and corrected (D + C). Hence, by comparing the RVE values computed with and without correction, we could verify if this approach could rectify the overall volume measurements.

Averaged real 3D DSC values and mean relative volume errors (RVE) are reported in Table 7 for each dataset and each class. The density plots of the predicted 2D DSC for the mid-slice from the original and corrupted datasets are displayed in figure 11. The density plots of the 150 predicted 3D DSC, for each dataset are displayed in figure 12 while the density plot of the RVEs for our three set- tings are displayed in figure 13. Regarding the impact of a poor quality segmentation in mid-slices on the 3D DSC values, we observed (Table 7) a substantial and significant 3D DSC reduction i.e. -18% for LVM (from 0.89 to 0.73), -12% for LVC, -27% for LVT and -35% for LVPM.

Changes regarding the relative volume errors were much larger (and all significant) i.e. +171% for LVM (from 0.07 to 0.19), +540% for LVC, +182% for LVT and +240% for LVPM. These results strongly supported that a relatively small 3D DSC change e.g. from 0.92 to 0.81 (0.7 has been proposed in the litterature as a good quality threshold for LVC) might correspond to a very large change in relative volume errors e.g. from 0.05 to 0.32 for LVC. The capacity of our 2D-based QC model to correctly detect the midslice errors was supported by the significant (according to Kolmogorov-Smirnov test) distributions shifts of both 2D DSC for each class (figure 11) and 3D DSC (figure 12). The effect was more pronounced for 2D DSC than for 3D DSC distributions. Table 7 andfigure

Discussion

In the present study, we reported a MISAQC tool designed to jointly predict 2D and 3D DSC values for the segmented structures in cardiac MR images. The tool performance was compared with a 3D-direct approach. Regardless of the segmented structure and the chosen backbone architecture and based on the accuracy and errors metrics, the computation of 3D DSC values based on 2D indices systematically outperformed the 3D-direct approach. The main contributions of our work are the following: • A mathematical solution allowing to get 3D DSC predictions from 2D DSC and mean volume similarity fraction (MVSF) predictions.

• A localization of segmentation errors at a 2D/slice level with a possibility to correct the effect of segmentation errors on clinical measurements.

• A significant improvement of results obtained from state of the art approaches (for 3D DSC predictions).

• A general approach, implementable for any MRI or computed tomography (CT) segmentation task, even if classes are unbalanced.

Comparison with existing works

In most of the previous studies which have addressed the issue of MISAQC with machine-learning-based methods, evaluations were carried on at the 3D level whereas the 2D-level was ignored [START_REF] Albà | Automatic initialization and quality control of large-scale cardiac mri segmentations[END_REF][START_REF] Audelan | Unsupervised quality control of image segmentation based on bayesian learning[END_REF][START_REF] Kohlberger | Evaluating segmentation error with-out ground truth[END_REF][START_REF] Robinson | Automated quality control in image segmentation: Application to the uk biobank cardiovascular magnetic resonance imaging study[END_REF][START_REF] Robinson | Real-time prediction of segmentation quality[END_REF][START_REF] Valindria | Reverse classification accuracy: Predicting segmentation performance in the absence of ground truth[END_REF]. Based on DSC computation at the 2D level, the present method yielded correlation factors higher than 0.965 regardless of the segmented structures. In a previous study conducted in lung images, the correlation between the predicted and the real DSC values was lower (r = 0.85) [START_REF] Kohlberger | Evaluating segmentation error with-out ground truth[END_REF], while values for liver and brain tumors were much lower (r = 0.54 and r = 0.69, respectively) [START_REF] Audelan | Unsupervised quality control of image segmentation based on bayesian learning[END_REF].

In the present study, the segmentation quality assessment was based on the DSC values, DSC being a continuous variable between 0 and 1. The corresponding evaluation allowed to score the segmentations quality on a continuous scale. Using 3D DSC quality thresholds, the reported accuracies were systematically superior to 96,5% for both the well represented (LVM and LVC) and the less represented classes (LVT and LVPM). As a matter of comparison, using a random forest classifier method for cardiac images, [START_REF] Albà | Automatic initialization and quality control of large-scale cardiac mri segmentations[END_REF] used a 4 mm error threshold and reported a 96% accuracy for LVM and LVC, two classes that are very well represented in cardiac images. The utilization of such a threshold for classes such as LVT, is highly questionable given their small and scattered structure. In addition, the 4 mm threshold enabled the distinction of only two groups of errors (lower and higher than 4 mm) which is a critical limitation for accuracy computations.

The MAE values we reported in the present study are similar to those reported by [START_REF] Robinson | Real-time prediction of segmentation quality[END_REF] who used a 3D-direct approach and a much larger dataset. More particularly, the MAE reported by these authors for LVM was 0.055 ± 0.064 as compared to 0.037 ± 0.035 in the present study using the QC ResNet 3D. Similarly, a value of 0.029 ± 0.039 was found in the present study for the LVC class using the QC ResNet 3D while Robinson et al reported a slightly larger value, i.e., 0.038 ± 0.040.

On that basis, these results strongly support the similarity between the two 3D-direct approaches. In addition, our 2D R approach clearly outperformed the 3D direct approach thereby indicating that a better performance can be achieved when MISAQC is initially performed at the 2D level. This new method was also very efficient in terms of computing time. The mean processing time for a subject was around 300 ms which permits real-time applications, on the contrary to other methods such as reverse classification accuracy where the corresponding reported time was 11 min (6000-fold larger) [START_REF] Robinson | Automated quality control in image segmentation: Application to the uk biobank cardiovascular magnetic resonance imaging study[END_REF]. We also addressed a limitation of previous works who used artificially degraded segmentations for learning and evaluation of their QC tools, and never confronted the automatic evaluation to human evaluation. For that we considered a clinical practice application using a real dataset having 1016 segmented subjects from the UK BioBank whose segmentation quality was assessed by a medical expert. Quality scores from our model were in very good agreement with those given by human operator, further supporting not only the robustness of our method but also it's ability to generalize to unseen data.

Advantages of a 2D-based multi-level MISAQC method

In the present method, training and inference were 2D-based whereas the quality assessment was performed at both 2D and 3D levels. In that respect, this MISAQC approach differs from those in which training, inference and quality evaluation have been performed at a 3D-level and provides multiple advantages which are listed below.

Possibilities opened by a multi-level QC

As indicated in section 5, the MISAQC performed at a slice-level allowed the detection of erroneous segmented slices (Table 2, figure 9 and 11). In addition, once detected, the corresponding erroneous measurement could be replaced by measurements from the neighbouring slices. In that respect volume errors could be drastically reduced (Table 7, figure 13). As a matter of example, the LVM mean RVE was 32% before the correction and was reduced by a factor five (5%) after the correction. This corrective process would not have been possible using a 3D-only quality evaluation. Being able to correct in real-time automatic clinical measurements can have an obvious diagnostic impact. Of interest, the impact of the mid-slice errors was moderate on the 3D DSC values and much larger on the volumetric errors (Table 7). This is a very important result indicating that the 3D DSC metric might not be sufficient in evaluating the quality of medical image segmentations and that morphological and/or functional metrics should be also considered. The additional 2D DSC evaluation at slice-level represents a clear advantage in that matter. Another application of the multi-level quality control would be to allow the clinicians to be able, in their diagnosis, to give more credit to the segmentation-derived measurements by being capable of directly visualizing the location of the segmentation errors and eventually correct them immediately. Using the feedback brought by our models, the clinician would also identify more easily and reliably the profile of data for which the model's predictions should be carefully checked.

Enhanced predictive accuracy

We showed that by reconstructing the 3D DSC from 2D indices one could obtain 3D-level and 2D-level quality evaluation of segmentations with a single method. This would not be possible in the case of a 3D-direct approach. This methodology also encourages the robustness of the 3D-level evaluation: since (as illustrated in Table 2) the 2D R versions can predict 2D indices with a very high predictive accuracy, the likelihood of large 3D DSC errors is substantially lower than what could be achieved using a 3D-direct approach.

Training and data efficiency

Annotation cost and computational cost are two recognized issues in the field of CNN-based segmentation. In that respect, using a 2D approach for a given dataset, the number of training examples is substantially increased [START_REF] Baumgartner | An exploration of 2d and 3d deep learning techniques for cardiac mr image segmentation[END_REF], without any additional annotation cost. Indeed, when working in a 3D framework, the patients segmentations' quality is expressed by a unique score for each class, whereas this information is enriched in the 2D context and multiplied by the number of slices. Having a higher variability in the training sample is highly beneficial for the convergence of the stochastic gradient algorithm. As a matter of example, in the promising study of [START_REF] Robinson | Real-time prediction of segmentation quality[END_REF] the training dataset was the UK Biobank database with thousand of available segmented 3D volumes. Such a gigantic database can be considered as exceptional and not commonly available. On that basis, a 2D dataset of an hundred segmented MRIs could be comparable to the UK BioBank database from the 3D point of view in terms of number of training examples. Obviously, computational costs are substantially reduced as well with our approach. This could be of particular relevance for computed tomography (CT) datasets for which the number of slices is of several hundreds: 3D approaches might face computational limits, whereas our approach will be easily implemented.

Quality control for small anatomical regions

Of particular interest, in addition to providing a quality control at both 2D and 3D levels, our method was insensitive to DSC values and size of the segmented structures. In other words, we were able to predict with the same accuracy both high and low DSC values. The issue of size structures was also addressed very efficiently. An interesting and challenging issue in the field of medical images segmentation is related to MISAQC for classes highly and poorly represented. Previous studies have mainly reported methods dedicated to the MISAQC of segmentations for relatively large and connected classes, such as the myocardium and blood cavity [START_REF] Audelan | Unsupervised quality control of image segmentation based on bayesian learning[END_REF][START_REF] Robinson | Automated quality control in image segmentation: Application to the uk biobank cardiovascular magnetic resonance imaging study[END_REF][START_REF] Robinson | Real-time prediction of segmentation quality[END_REF]. The suitability and performance of these methods for small and scattered anatomical structures such as trabeculations or papillary muscles (illustrated in figure 1), have not been addressed. The present results clearly demonstrate that our 2D R method is suitable for both large and small classes. The corresponding MAEs were similar (QC ResNet 2D R : 0.016 (LVM), 0.011 (LVC), 0.022 (LVT) and 0.042 (LVPM)), illustrating that the CNNs were able to learn and extract relevant features for both types of classes. Our method was less efficient for the LVPM class and so most likely because papillary muscles do not appear in all slices.

Limitations

A few limitations should be acknowledged in the present work. One could wonder about the realistic features of the synthetic database we generated. The methodology we used to generate this database was comparable to what has been previously reported [START_REF] Robinson | Real-time prediction of segmentation quality[END_REF]. Actually, we replaced the random forests approach by CNNs and this slight change is not expected to introduce major changes. More importantly, the external validation of our model on the UK BioBank dataset supports the realistic nature of this synthetic dataset as the model was able to generalize on a new data source. We used a single type of network to build the MISAQC dataset and one could wonder whether other networks could have produced other "types" of segmentation errors. We think that the key issue for this dataset construction is to generate a sufficient number of segmentations fairly distributed over the entire range of possible DSCs (from 0 to 1).

The performance of our approach was assessed using a rather small dataset. However, as it relies on a 2D analysis, it should be less sensitive to the dataset size than a 3D-based deep learning approach. This expectation is supported by the promising comparison between our results and those from [START_REF] Robinson | Real-time prediction of segmentation quality[END_REF] with a smaller dataset, we reported lower MAE values.

Conclusion

In the present study, a multi-level (2D and 3D) deep-learning-based real-time automated quality control method for cardiovascular MR image segmentations was designed and the corresponding performance was assessed. Based on the MAE values and classification accuracy, it was clearly demonstrated that the proposed method was equally efficient for large and small cardiac anatomical structures. The QC provided by our method was shown to be consistent with scores generated by trained cardiologists. The 2D-based structure can be trained with a modestly sized dataset while enabling a very accurate real-time automatic quality prediction at both 2D and 3D levels. We highlighted its' possible applications such as successfully rectifying erroneous clinical measurements derived from medical image segmentations.

Fig. 1 .

 1 Fig. 1. First row: a CMR image from our testset and related manual segmentation; the four left ventricle structures are depicted illustrating the difference in size and shape for the left ventricle trabeculations (LVT) and left venctricle papillary muscles (LVPM) structures compared to left ventricle myocardium (LVM) and left venctricle cavity (LVC). Second row: a related automated segmentation to evaluate; the displayed scores represent for each class the real DSC on the left and the predicted DSC with our method on the right.

Fig. 2 .

 2 Fig. 2. Illustration of MISAQC dataset construction. First, 4290 cardiac MR images were segmented by an expert and used to train segmentation models in order to produce segmentations of varying quality of the 4290 MR images. Then, the DSC values were computed for those segmentations based on the ground truth manual segmentations. Finally, LVT DSC-based subsampling was applied before splitting the data into training/validation and test sets.

Fig. 3 .

 3 Fig. 3. Illustration of the methodology at inference time. For a given segmented MRI to evaluate, the inference is done slice-wise by sequentially providing the segmentation and related image as the input of the two models. The predicted 2D DSCs are stored and provide a quality evaluation for all slices and classes. When slice-wise inference is over, the stored predictions are combined to compute the 3D DSC prediction for all classes. Errors can be detected at every scale.

Fig. 4 .

 4 Fig. 4. Display of the 4 used architectures, for each backbone type the numbers in red are only present for the 3D version. Downsampling was always performed with a max pooling layer while the upsampling (only in QC U-Net) was performed with transposed convolutions. The pooling size can be deduced from the difference in dimensions between consecutive layers. Transposed convolution kernels were always 2x2 (2x2x1 in 3D version). QC ResNet 2D and 3D had 2 and 3 million parameters respectively, while QC Unet 2D and 3D had 19 and 21 million parameters respectively.

Fig. 6 .

 6 Fig. 6. Boxplots of the 3D DSC absolute errors (AE) for each class and method.

Fig. 7 .

 7 Fig. 7. Scatter plots of the true 3D DSC values against the predicted values for all classes except BG. The R scores and red lines are derived from a fitted regression model. The model used here is QC U-Net 2D R .

Fig. 8 .

 8 Fig. 8. 3D DSC MAE (averaged over all LV classes) as a function of the percentage of total training data used during training. The model used here is QC U-Net 2D R .

Fig. 9 .

 9 Fig. 9. Density plots of the predicted 2D DSC for the QC 1 and QC 2 and QC 3 groups in clinical evaluation. Within each class, all distributions were significantly different except QC 2 vs QC 3 for LVT and LVPM.

Fig. 11 .

 11 Fig. 11. Density plots of the predicted 2D DSC for the mid-slice segmentation, before and after deterioration, on 150 subjects in synthetic data. Within each class, all distributions were significantly different.

  13 indicated that the shift in measurement errors distributions induced by the mid-slice errors could be almost totally suppressed using our proposed correction. This was further confirmed by Kolmogorov-Smirnov test as the RVE w/o D and RVE D distributions were always found significantly different while the RVE w/o D and RVE D + C distributions were not (p-values: 0.99 for LVM, 0.95 for LVC, 0.89 for LVT and 0.72 for LVPM).

Fig. 12 .

 12 Fig. 12. Density plots of the predicted 3D DSC for 150 3D segmentations, before and after deterioration in synthetic data. Within each class, all distributions were significantly different.

Table 1 . Baseline Patient Characteristics for the initial dataset.

 1

		Characteristics
	Parameter	HCM	DCM	ETCM	Healthy	Overall
	Total no.	92	39	59	259	449
	Age (y)	53 ± 15	54 ± 12	44 ± 16	40 ± 17	45 ± 17
	No. men	59	21	37	147	264
	LVEDV (mL)	151.5 ± 38.1 226.5 ± 86.7 189.0 ± 90.3 146.4 ± 42.3 160.0 ± 60.3
	LVEDV-to-BSA (mL/m )	80.0 ± 17.4 118.8 ± 41.0 101.4 ± 39.0 80.5 ± 19.7	86.5 ± 27.8
	LVMM (g)	169.8 ± 53.8 163.9 ± 48.1 124.3 ± 66.8 113.9 ± 34.3 131.0 ± 51.2
	LVMM-to-BSA (g/m )	89.6 ± 26.2	86.2 ± 22.4	66.6 ± 29.2	62.3 ± 14.7	70.5 ± 23.5
	Trabeculation mass (g)	11.3 ± 5.8	17.0 ± 10.0	18.9 ± 13.3	9.5 ± 4.7	11.7 ± 7.9
	Trabeculation mass-to-BSA (g/m )	6.0 ± 3.0	8.8 ± 4.7	10.1 ± 6.3	5.2 ± 2.4	6.3 ± 3.9
	Trabeculation mass-to-TMM (%)	6.0 ± 2.4	8.6 ± 3.8	12.7 ± 5.2	7.1 ± 2.6	7.6 ± 3.7
	0.4], [0.4, 0.5], [0.5, 0.6], and [0.7, 1]. The number of
	samples taken from each interval is the minimum number
	observed within these intervals. Hence, 708 3D segmen-
	tations were randomly selected in each bin, resulting in a
	total of 4248 3D segmentations, corresponding to 40 363
	2D segmentations. The final dataset was randomly split
	into training/validation (90%) and testing (10%) samples,
	stratifying over the LVT 3D DSC for both.

Table 2 . DSC MAEs, accuracies and correlation coefficients at slice-level for each 2D R version.

 2

				MAE
	Method	LVM		LVC		LVT	LVPM
	QC ResNet 2D R	0.025 ± 0.027	0.021 ± 0.025	0.047 ± 0.048	0.023 ± 0.044
	QC U-Net 2D R	0.023 ± 0.024	0.019 ± 0.022	0.047 ± 0.051	0.020 ± 0.040
	p-values	3.62e-05	2.18e-07	0.47		6.01e-30
			Accuracy and correlation
		LVM		LVC		LVT	LVPM
	Method	Acc.	Corr.	Acc.	Corr.	Acc.	Corr.	Acc.	Corr.
	QC ResNet 2D R	96.1 %	0.991	97.3%	0.992 92.8% 0.960	98.5%	0.991
	QC U-Net 2D R	96.6 % 0.993 97.4 % 0.994 92.7 % 0.957 98.8 %	0.993

Table 3 . DSC MAEs at subject-level for each backbone and method. Regardless of the backbone the 2D R method consis- tently outperforms the equivalent 3D direct approach.

 3

			QC ResNet backbone
	Method	BG	LVM	LVC	LVT	LVPM
	3D	0.011 ± 0.022	0.037 ± 0.035	0.029 ± 0.039	0.044 ± 0.038	0.055 ± 0.049
	2D p-values	9.96e-30	2.47e-34	3.02e-31	5.55e-25	5.12e-07
			QC U-Net backbone
	Method	BG	LVM	LVC	LVT	LVPM
	3D	0.010 ± 0.018	0.036 ± 0.035	0.026 ± 0.031	0.045 ± 0.041	0.042 ± 0.039
	2D p-values	2.42e-24	6.67e-34	7.57e-26	2.08e-21	0.02

R 0.004 ± 0.005 0.016 ± 0.028 0.011 ± 0.016 0.022 ± 0.022 0.042 ± 0.053 R 0.004 ± 0.006 0.016 ± 0.028 0.012 ± 0.017 0.023 ± 0.023 0.040 ± 0.054

Table 4 . Accuracies and correlation coefficients at subject-level for backbone and method.

 4

				QC ResNet backbone
		LVM		LVC		LVT		LVPM
	Method	Acc.	Corr.	Acc.	Corr.	Acc.	Corr.	Acc.	Corr.
	3D	92.7 %	0.983	94.7 %	0.977	90.7 %	0.955 96.5 % 0.961
	2D R	97.7% 0.992 97.0% 0.996 97.2% 0.987 97.2% 0.965
				QC U-Net backbone
		LVM		LVC		LVT		LVPM
	Method	Acc.	Corr.	Acc.	Corr.	Acc.	Corr.	Acc.	Corr.
	3D	94.4 %	0.982	96.0 %	0.986	89.7 %	0.951 97.7 % 0.971
	2D R	97.5 % 0.992 97.0 % 0.996 96.5 % 0.987 97.7 % 0.965

Table 6 . Mean processing time for each network at subject-level and at slice-level additionally for 2D R versions.

 6

	Inference time
	Method	Subject-wise Slice-wise
	QC ResNet 3D	100 ms	∅
	QC U-Net 3D	120 ms	∅
	QC ResNet 2D R	314 ms	32 ms
	QC U-Net 2D R	312 ms	31 ms
	0.042 (0.053) for LVPM.

Table 7 . Additional analyses: cases with and without deteriorated mid-slice segmentation.

 7 ± 0.07 0.07 ± 0.05 0.19 ± 0.14 0.07 ± 0.05 LVC 0.92 ± 0.02 0.81 ± 0.09 0.05 ± 0.05 0.32 ± 0.29 0.05 ± 0.05 LVT 0.62 ± 0.08 0.45 ± 0.11 0.23 ± 0.16 0.65 ± 1.06 0.23 ± 0.15 LVPM 0.79 ± 0.11 0.51 ± 0.12 0.15 ± 0.13 0.51 ± 0.40 0.16 ± 0.14

			3D DSC and RVE
	Class	3D DSC w/o D 3D DSC D RVE w/o D	RVE D	RVE D + C
	LVM	0.89 ± 0.26	0.73

Acknowledgments

Steffen E. Petersen (SEP) acknowledges support from the National Institute for Health Research (NIHR) Biomedical Research Centre at Barts. SEP acknowledge the British Heart Foundation for funding the manual analysis to create a cardiovascular magnetic resonance imaging reference standard for the UK Biobank imaging resource in 5000 CMR scans (www.bhf.org.uk; PG/14/89/31194). SEP acknowledges support from the "SmartHeart" EP-SRC programme grant (www.nihr.ac.uk; EP/P001009/1). SEP and ER also acknowledge support by the London Medical Imaging and Artificial Intelligence Centre for Value Based Healthcare (AI4VBH), which is funded from the Data to Early Diagnosis and Precision Medicine strand of the government's Industrial Strategy Challenge Fund, managed and delivered by Innovate UK on behalf of UK Research and Innovation (UKRI). Views expressed are those of the authors and not necessarily those of the AI4VBH Consortium members, the NHS, Innovate UK, or UKRI. SEP provides consultancy to and owns stock of Cardiovascular Imaging Inc, Calgary, Alberta, Canada. Other authors did not receive any financial or material support from any industrial company.