Benoit Goutorbe
email: benoit.goutorbe@inrae.fr

Anne-Laure Abraham

Ghislain Bidaut

Philippe Halfon

Mahendra Mariadassou

Anne Plauzolles

Sophie Schbath

Characterizing the limits of shallow shotgun metagenomics for taxonomic pro ling of human gut microbiota in clinical studies

Keywords: microbiota, metagenomics, shotgun metagenomics, shallow shotgun metagenomics microbiota, metagenomics, shotgun metagenomics, shallow shotgun metagenomics

Background: Shallow shotgun metagenomics (SSM) has been recently suggested as a promising strategy to study human microbiota, providing nearly identical taxonomic profiles to deep shotgun metagenomics but at a sequencing cost as low as that of metabarcoding. To help clinical researchers determine whether shallow sequencing is appropriate for their projects, it is crucial to ascertain the accuracy of the information it provides, compared to deep sequencing. Here, we design a mapping-based workflow to build taxonomic profiles from SSM data and assess its accuracy at varying sequencing depths at both sample and cohort levels using extensive simulations and several public data sets.

Results: To identify genuinely present species and spuriously identified ones, we propose a novel data-driven filtering method based on machine learning techniques that largely outperforms basic filtering strategies based on predefined thresholds, resulting in reliable taxonomic profiles at different sequencing depths, ranging from 50K to 10M reads/samples. Up to 90% of species with relative abundances higher than 4.10 -4 were recovered correctly at 500K reads/sample showing that only information about rare taxa is lost at shallow depths. Furthermore, our results clearly show that SSM is able to correctly recover relevant biological signal from the confidently identified taxa, such as differences between groups of patients and diagnosis-like classification.

Conclusions: This study confirms that SSM is suitable for clinical research on human gut microbiota. We recommend that researchers should consider moving from 16S to SSM to limit biases in taxonomic profiles, or moving from deep to shallow sequencing, when functional analyses are not the main focus, to reduce costs and be able to include more patients in research projects.

Background

High throughput sequencing (HTS) unraveled the quantity, diversity and complexity of host-microbiota interactions in health and human diseases by allowing culturefree analyses of microbial ecosystems [START_REF] Fan | Gut microbiota in human metabolic health and disease[END_REF] [2] [START_REF] Clemente | The role of the gut microbiome in systemic inflammatory disease[END_REF]. Metabarcoding, which consists in targeted sequencing of a phylogenetic marker (often the 16S rRNA gene for bacteria), has been widely used, as it allows to assess the biodiversity within and between samples, and to obtain an approximate taxonomic identification (most frequently down to the genus level) of the microorganisms present in a sample [START_REF] Golob | Evaluating the accuracy of amplicon-based microbiome computational pipelines on simulated human gut microbial communities[END_REF]. Besides, whole genome sequencing (WGS, also referred as shotgun sequencing) allows deeper taxonomic resolution (down to species, or even strain level) [START_REF] Ellegaard | Beyond 16S rRNA Community Profiling: Intra-Species Diversity in the Gut Microbiota[END_REF] [START_REF] Olm | inStrain profiles population microdiversity from metagenomic data and sensitively detects shared microbial strains[END_REF], functional profiling (identification and quantification of genes, metabolic pathways) [START_REF] Abubucker | Metabolic Reconstruction for Metagenomic Data and Its Application to the Human Microbiome[END_REF], and de novo assembly of uncultured organism genomes as Metagenome-Assembled Genomes (MAGs) [START_REF] Parks | Recovery of nearly 8,000 metagenome-assembled genomes substantially expands the tree of life[END_REF]. Even though it is often claimed that human microbiota research would need to head toward WGS to deeply understand host-microbiota interaction [START_REF] Ellegaard | Beyond 16S rRNA Community Profiling: Intra-Species Diversity in the Gut Microbiota[END_REF] [START_REF] Yen | Metagenomics: a path to understanding the gut microbiome[END_REF] [10] [START_REF] Durazzi | Comparison between 16S rRNA and shotgun sequencing data for the taxonomic characterization of the gut microbiota[END_REF], metabarcoding is still widely used, notably because of its substantially lower sequencing and data processing costs. Indeed, researchers often favour the number of samples to include many patients and/or have multiple samples per patient (longitudinal approach), rather than the amount of information per sample.

Shallow shotgun metagenomics (SSM) have been recently suggested as an alternative middle route [START_REF] Hillmann | Evaluating the Information Content of Shallow Shotgun Metagenomics[END_REF], which is both cost-competitive with metabarcoding and as informative as deep WGS for taxonomic profiling. Whereas tens of millions of reads per sample are typically used to characterize human gut microbiota samples with standard (deep) shotgun sequencing [START_REF] Qin | A human gut microbial gene catalog established by metagenomic sequencing[END_REF] [START_REF] Quince | Shotgun metagenomics, from sampling to analysis[END_REF], SSM typically deals with fewer than 1M reads/sample, thus drastically reducing sequencing costs. Previous works suggest that species level taxonomic profiles obtained by mapping reads on reference genomes were highly similar to deep WGS, and that the resulting α-diversity metrics were barely impacted by limiting the sequencing depth to ∼ 500K -1M reads/sample [START_REF] Hillmann | Evaluating the Information Content of Shallow Shotgun Metagenomics[END_REF] [15] [START_REF] Gweon | The impact of sequencing depth on the inferred taxonomic composition and AMR gene content of metagenomic samples[END_REF] [START_REF] Cattonaro | Do you cov me? Effect of coverage reduction on metagenome shotgun sequencing studies[END_REF]. Furthermore, mapping on reference genomes was more efficient than on a marker genes catalog in the context of SSM [START_REF] Santiago-Rodriguez | Metagenomic Information Recovery from Human Stool Samples Is Influenced by Sequencing Depth and Profiling Method[END_REF]. By contrast, the depth required for functional analysis depends on the level of granularity of the analysis: 500K reads/sample may be sufficient to identify KEGG orthology groups [START_REF] Hillmann | Evaluating the Information Content of Shallow Shotgun Metagenomics[END_REF], but 3M to 5M reads/sample are necessary to accurately detect genes and pathways [START_REF] Santiago-Rodriguez | Metagenomic Information Recovery from Human Stool Samples Is Influenced by Sequencing Depth and Profiling Method[END_REF] [START_REF] Treiber | Pre-and post-sequencing recommendations for functional annotation of human fecal metagenomes[END_REF], and very deep sequencing up to ∼ 60 -80M reads/sample is needed to study antimicrobial resistance genes [START_REF] Gweon | The impact of sequencing depth on the inferred taxonomic composition and AMR gene content of metagenomic samples[END_REF] [START_REF] Zaheer | Impact of sequencing depth on the characterization of the microbiome and resistome[END_REF]. However, complementary investigations on the reliability of taxonomic profiles constructed by mapping SSM reads on a catalog of representative genomes are needed. In particular, as we deal with critically low sequencing depths, it is essential to retrieve as much information as possible from each read, including ambiguous reads that map to several genomes, and to identify genomes present in low abundances. To the best of our knowledge, previous works used catalogs that did not include MAGs (Metagenome-Assembled Genomes), limiting resulting profiles to cultured organisms. Huge metagenomic efforts have been made to study the human gut microbiota, allowing to build exhaustive catalogs, that gather both cultured and uncultured organisms, expanding taxonomic resolution for following metagenomic studies [START_REF] Almeida | A unified catalog of 204,938 reference genomes from the human gut microbiome[END_REF]. In addition, there is still limited knowledge about the potential loss of biological signal recovery in clinical metagenomics data sets when switching from deep WGS to SSM.

In the present study, we aim (i) to build a mapping-based workflow for SSM data to build accurate taxonomic profiles, that uses state of the art reference genome database, (ii) to evaluate the effect of sequencing depth on taxonomic profiles and, and (iii) underlying signal recovery (e.g. stratification of patients into groups) for a clinical study. To do so, we used extensive simulations and analyzed several publicly available data sets.

Methods

Data

Simulated data sets.

We retrieved taxonomic profiles of 100 human gut microbiotas from Qin 2014 [START_REF] Qin | Alterations of the human gut microbiome in liver cirrhosis[END_REF] through curatedMetagenomicData [START_REF] Pasolli | Accessible, curated metagenomic data through ExperimentHub[END_REF]. The samples have a richness of 98 ± 15 species per sample, and species' relative abundances range from 5.10 -1 down to 10 -6 (geometric mean 5.10 -4 ). These profiles were given using the NCBI's taxonomy and were converted into UHGG's taxonomy [START_REF] Almeida | A unified catalog of 204,938 reference genomes from the human gut microbiome[END_REF] v1.0 by choosing the UHGG species with the taxonomic assignation closest to the NCBI one (if several NCBI species were tied, one was chosen randomly), resulting in profiles with the exact same complexity, approximately the same phylogenetic composition but including some uncultured organisms (MAGs). UHGG genomes are clustered into "species clusters" (thereafter referred to as species), that share at least 95% of identity on 30% of the genomes. One representative genome is chosen in each cluster for inclusion in the mapping catalog. We generated profiles using a randomly selected genome of each targeted species, in order to mimic the case where the strain we observe is not the one present in the catalog. For each of the 100 samples, we simulated 10M paired end reads, using Grinder [23] v0.5.3 (length of 2*125bp, insert size normally distributed with an average of 500bp and standard deviation of 50 bp without sequencing error) and subsampled the read sets at 5 M , 1 M , 500 K, 100 K, 50 K and 10 K reads/sample.

Real data sets.

We used data from 3 clinical studies, for a total of N = 439 samples covering patients from several continents and clinical conditions (healthy patients, hepatic diseases at different stages, cancer patients). Loomba et al. (2017) [START_REF] Loomba | Gut Microbiome-Based Metagenomic Signature for Non-invasive Detection of Advanced Fibrosis in Human Nonalcoholic Fatty Liver Disease[END_REF] compared the gut microbiota of 86 patients suffering from hepatic diseases at different stages (N = 14 fibrosis vs N = 72 NAFLD). [START_REF] Matson | The commensal microbiome is associated with anti-PD-1 efficacy in metastatic melanoma patients[END_REF] [START_REF] Matson | The commensal microbiome is associated with anti-PD-1 efficacy in metastatic melanoma patients[END_REF] compared, among 39 patients having metastatic melanomas, those who responded to anti-PD-1 immunotherapy (N = 15) and those who did not (N = 24). [START_REF] Qin | Alterations of the human gut microbiome in liver cirrhosis[END_REF] [START_REF] Qin | Alterations of the human gut microbiome in liver cirrhosis[END_REF] compared patients having liver cirrhosis (N = 169) and a group of healthy controls (N = 145), with a discovery and a validation cohort for both groups. We analyzed these data sets at full depth and subsampled them to mimic shallow sequencing in the remainder, as described above.

Bioinformatics pipeline

Reads were pre-processed using trimmomatic [START_REF] Bolger | Trimmomatic: a flexible trimmer for Illumina sequence data[END_REF] v0.39 to remove low quality reads (with average quality below 30) and reads shorter than 80 nucleotides. For real data sets, reads were also mapped to the human genome (hg38) to filter out host contamination. Remaining reads were then mapped to the UHGG catalog [START_REF] Almeida | A unified catalog of 204,938 reference genomes from the human gut microbiome[END_REF] v1.0 using bwa mem [START_REF] Li | Aligning sequence reads, clone sequences and assembly contigs with BWA-MEM[END_REF] v0.7.17 local aligner, with option -h50 to allow up to 50 reported hits. Mapping files were then processed using samtools [START_REF] Li | The Sequence Alignment/Map format and SAMtools[END_REF] v1.9 and custom scripts developed under Python v2.7.13.

Ambiguous reads (i.e. reads that map to several genomes) occur frequently when mapping reads to a catalog of reference genomes (26% and 42% of the mapped reads in simulated and real data sets respectively), due to highly conserved genes and mobile elements notably. Thus, we split mapped reads into unambiguous reads (i.e. mapping to one genome only), and ambiguous reads. For each genome, we computed the reads count (RC) and the fraction of the genome covered (FC) by at least one read, using either all reads or unambiguous reads only (uRC and uFC), as well as a specificity ratio (SR) defined by the number of unambiguous reads divided by the total number of reads mapped to this genome (uRC/RC).

In order to estimate the species' s relative abundances, we first compute the representative genomes' average coverage C s = 1 ℓs i r i,s , with ℓ s being the length of the representative genome of species s and r i,s the length of read i that is unambiguously mapped to s. We then obtain the relative abundance by normalizing across species to sum up to 1: A s = Cs j Cj . We refine this estimation by reallocating the ambiguous reads, randomly assigning them to one of their hits with a probability proportional to A s .

Simulations analysis

Direct mapping of short reads on reference genomes produces many false positives (genomes covered by reads but not present in the sample) that need to be filtered out. We used simulated profiles, with known composition, to determine the most efficient way to classify the genomes into true positives (TP) and false positives (FP). In order to assess methods and compare them to each other, we computed the area under the receiver operating characteristic (ROC) curve (AUC) for this classification task, using evabic R package (https://github.com/abichat/evabic). We also implemented an automated threshold search, that limits the false discovery rate (F DR = F P T P +F P ) to a maximum of 0.1, and compared false negative (FN) rates at this threshold across methods and sequencing depths.

We first evaluated how genomes' features (RC, uRC, FC, uFC and SR) can be used independently to classify the genomes, and then combined them to train classifiers. We used logistic regression (LR), linear discriminant analysis (LDA) and random forests (RF), to perform classification, with uRC, uFC, SR and total sequencing depth as input features. Finally, we used a 4-fold cross validation procedure to evaluate the performance of these methods and determine suitable thresholds for each method and sequencing depths.

Real data sets analysis

We analyzed real data sets using (1) RF-based filters fitted on the simulations data and thresholds that control FDR at each sequencing depth, and (2) a basic filtering that discards all species with a relative abundance beyond 10 -4 , FC beyond 10 -2 or uFC beyond 10 -4 . This filtering is inspired by the one used by Santiago-Rodriguez et al. [START_REF] Santiago-Rodriguez | Metagenomic Information Recovery from Human Stool Samples Is Influenced by Sequencing Depth and Profiling Method[END_REF]. Note that it is a quite permissive threshold due to the low sequencing depths to which it is applied.

We evaluated α-diversity using species richness and Shannon diversity, and βdiversity using Jaccard distance and Bray-Curtis dissimilarity index, computed with the R package phyloseq [START_REF] Mcmurdie | phyloseq: An R Package for Reproducible Interactive Analysis and Graphics of Microbiome Census Data[END_REF]. In order to assess the impact of sequencing depth on taxonomic profiles, we evaluated the correlation between subsampled and deep αdiversity measures using Spearman and Pearson correlation as well as the correlation between species' relatives abundances at full and shallower depths. Additionally, we measured the distance between low depth samples and their full depth counterparts. Finally, for each data set, we evaluated whether differences between groups of interest were preserved at lower sequencing depths using the following criteria :

• difference in α-diversity between groups, through a Wilcoxon test on the aforementioned metrics,

• structure in the β-diversity matrix, through a PERMANOVA analysis of the aftermentioned dissimilarity indices,

• biomarker discovery, using a Wilcoxon test with Benjamini-Hochberg correction to recover differentially abundant species,

• patients' classification in their groups of interest, using random forests trained on taxonomic profiles together with species richness and Shannon's diversity, and performing an iterative feature selection step as described by Loomba et al. [START_REF] Loomba | Gut Microbiome-Based Metagenomic Signature for Non-invasive Detection of Advanced Fibrosis in Human Nonalcoholic Fatty Liver Disease[END_REF], with 10 repetitions at each step to take into account the variability of the classification method. For Loomba-2017, we asked authors for patients age and BMI, and imputed missing values with the mean value of the cohort. For Qin-2014, we performed one classification with all patients, and another classification where the discovery cohort was used for training and the validation cohort for testing. In order to perform unbiased comparisons of p-values and AUCs across sequencing depths, we used only samples having at least 10M high quality reads per sample for Loomba-2017 (N = 77) and Matson-2018 (N = 39) and only 5M reads for Qin-2014 (N = 236) as using the same 10M depth would have reduced the cohort analyzed down to 172 patients.

Results

A tailored filtering strategy greatly improves species recovery

We first used the 100 simulated metagenomes to design an appropriate pipeline to build taxonomic profiles from SSM reads mapped on a reference genomes catalog. Our raw mapping data from simulations shows that a small number of reads (8% of unambiguously mapped reads) are mapped to a large number of unexpected genomes, not part of the original profile used for read simulation. Overall, those numerous unexpected genomes have a low coverage but result in a great number of false positives if no filtering is applied (FDR = 0.92). The basic filtering approach, which discards the rarest and least covered species (see above for details), yields an overall FNR = 0.44 and FDR = 0.46. The different features available have contrasted discriminatory powers as measured by AUC: 0.76 for read counts (RC) and 0.87 for fraction covered (FC). The AUC increases when using only unambiguous reads: 0.84 for uRC and 0.90 for uFC (see Table 1 for details). We therefore used the latter two features to build optimized filters. As seen on Figure 1A, a thresholding strategy based on uFC and/or uRC values (corresponding to horizontal and/or vertical lines) to discriminate TPs and FPs, would be suboptimal as it would miss the long tail of genomes with low uRC but comparatively high uFC values. These observations motivated the development of a data-driven classifier that could leverage this pattern. We trained several classifiers (LDA, Logistic regression, Random forest) using uRC, uFC, SR, as well as sequencing depth as predictors to predict genomes' status (present or absent).

We can see in Table 1 that data-driven classifiers largely outperform basic filtering. LDA and LR perform similarly, and yield nearly identical results in training and testing samples in the cross validation process, highlighting very good generalization capabilities. RF appears to be the best method, yielding a nearly perfect classification in training sets, and still outperforming others in the testing sets. RF was thus used in the rest of the work. Table 1 also shows that when choosing a threshold that limits the FDR at 0.1, data-driven classifiers have a drastically lower FN rates (∼ 0.3 -0.4) than the ones obtained with independent thresholding on each genomes features (∼ 0.6 -0.9).

The FN values are still quite high and may correspond to a subset of species that are intrinsically difficult to identify. To further investigate this issue and characterize the information loss induced by SSM, we plotted the distribution of simulated relative abundances of species that were absent in profiles with respect to the sequencing depth, as seen on Figure 1B, using RF-based filtering. We can clearly see the inflation of FN when lowering sequencing depth, but we can also notice that, as expected, the populations that are lost are relatively rare. For instance, at 500K reads/sample, 90% of species with an abundance greater than 4.10 -4 were detected. In comparison, 90% of the species with an abundance greater than 2.10 -4 were detected at 1M reads/sample, and this value was down to 3.10 -5 at 5M reads/sample.

While focusing on TPs, we observed that Pearson correlation between simulated and estimated species relative abundances increases from ρ = 0.81 to ρ = 0.91 if we add a step of reallocation of ambiguous reads. Therefore, we used profiles after reallocation of the ambiguous reads for the rest of the work.

Shallow shotgun metagenomics accurately reconstruct taxonomic profiles for abundant taxa

Applying the RF-based filters on real data sets results in profiles with an average diversity of 128 ± 66 species per sample at full depth. The diversity gradually decreases as the sequencing depth decreases, shrinking down to 45 ± 21 at 500K reads/sample for example (fig 2A). In contrast to the observed species count, the Shannon diversity index (fig 2B) is much less impacted by sequencing depth, indicating that the species lost at low sequencing depths are mostly rare ones, as expected. Down to 500K reads per sample, the correlation between full depth and subsampled Shannon indices remains very high. Distances between subsamples and their reference, defined as the corresponding sample at full depth, gradually increase when decreasing the sequencing depth (fig 2C), and show high replicability across data sets.

In comparison, basic filtering is more permissive, producing profiles with increased diversities and less affected by reduced sequencing depths (fig 2D ,E). Distances to the reference are also smaller than with RF-based filters (fig 2C ,F).

These results show that taxonomic profiles are not impacted by the sequencing depth for most abundant species, regardless of the filtering strategy used. For the rarest species, we can see that reducing sequencing depth gradually decreases the number of detected species, especially when using a filtering strategy that ensures few spuriously detected species.

Shallow shotgun metagenomics recover biological signal

Here, we assess the impact of the sequencing depth on biological signal recovery, i.e. stratification of patients into groups, in three different data sets covering several diseases and geographical origins. In each study, we track differences in the microbiota composition between two groups of patients according to a clinical condition (see Methods section for details). We aim to see whether the results obtained in deep shotgun sequencing are preserved when switching to SSM.

As expected from previous results, differences in α-diversity between groups are maintained using shallow sequencing: p-values were concordant across sequencing depths (see Fig. 3A), with a strong difference between groups in Qin-2014 data set at all depths, a slight but not significant difference between groups in Loomba-2017 and no difference between groups in Matson-2018. PERMANOVA analysis leads to similar results (Fig. 3B), showing that the structure of the distance matrix between samples is only marginally impacted by reducing sequencing depths. The p-value between healthy and sick patients in Loomba-2017 slightly increases when reducing the sequencing depth down to 500K but drastically increases when moving to 100K reads/sample and below, revealing that some key species for the stratification of patients may be lost at very low sequencing depths.

Moreover biomarkers discovery (Supplementary File 1) shows that the number of differentially abundant taxa identified in the data sets strongly depends on the sequencing depth, with abundant biomarkers easier to recover at shallow depths than rare ones. Indeed, out of the 6 biomarkers identified in Loomba-2017 with deep sequencing, the 4 markers with a geometrical mean abundance greater than 10 -2 were recovered at 500K reads/sample. Similarly, in Qin-2014, out of the 56 biomarkers identified with deep sequencing , the 2 markers with a geometrical mean abundance greater than 10 -2 and 12 of 37 markers with a geometrical mean abundance between 10 -2 and 10 -3 were recovered at 500K reads/sample. None of the rares ones, with a geometric mean abundance lesser than 10 -3 was recovered in both datasets. (see supplementary file 2 for details).

Finally, we performed classification of patients using RF on the taxonomic profiles in Loomba-2017 and Qin-2014 (see Fig. 3C). On Loomba-2017, we observed that, under 5M reads/sample, AUC gradually decreases as sequencing depth becomes smaller, due to the loss of some taxa that are crucial for the classification. On Qin-2014 data set, we could perform a very good classification on both discovery and validation cohorts even at low sequencing depths, with performance gradually decreasing under 5M reads/sample.

Results obtained by basic filtering lead to the same conclusion regarding the impact of the sequencing depth on patients stratification, as shown in supplementary file 3. In comparison with RF-based filtering, p-values tend to be lower and less impacted by sequencing depth on Loomba-2017, and AUCs upper. Basic filters are more permissive, thus more species are identified, especially rare ones at low sequencing depths. As seen in the simulations, they also include many spuriously identified species, but it seems that patients stratification is not widely impacted by the noise created by those spurious species.

Discussion

Our simulations highlighted the need for an efficient filtering strategy while mapping shallow shotgun metagenomic reads against a reference catalog to reconstruct taxonomic profiles. In order to have reliable results for every sequencing depth, depth-dependent thresholds were applied. This step is crucial to prevent misleading interpretations and to provide a trustworthy biological knowledge. Controlling the false discovery rate (FDR) in taxonomic profiles had the direct consequence of decreasing the number of identified species, especially at low sequencing depths. We observed a great benefit of random forest-based filters, and to a lesser extent of other machine learning-based models tested, in comparison with simple filters based on species features independently (read counts and fraction covered, considering all reads and unambiguous reads only), allowing to identify more species and rarer ones for equivalent FDR. On the three real data sets considered, our analysis showed that differences between groups of patients observed at full depth were still recovered at low sequencing depths. Permissive and depth-independent filtering, as performed in previous works on SSM, allowed a little improvement in structure recovery compared to our stringent random forest-based filters. Indeed, these structures were less sensitive to the noise introduced by spuriously identified species in the profiles using basic filtering, than to the removal of key species induced by our stringent random forest-based filter. However, for trustfulness and interpretability of results, downstream analysis should be based on reliably identified taxa even if it induces a moderate loss of biological signal. Regardless of the filtering strategy, our analysis shows that SSM is sufficient to recover most of the differences between groups of patients discovered with deep sequencing, provided that the differences are not based on rare species.

Our simulations led us to develop stringent filters, especially at low sequencing depths, resulting in profiles with limited species richness. Other studies reported a smaller effect on species richness while using similar mapping strategies but less stringent filters [START_REF] Hillmann | Evaluating the Information Content of Shallow Shotgun Metagenomics[END_REF] [START_REF] Santiago-Rodriguez | Metagenomic Information Recovery from Human Stool Samples Is Influenced by Sequencing Depth and Profiling Method[END_REF]. Unlike other published papers, the catalog we used contained metagenome-based assembled genomes, which are often incomplete and have contamination. It enables the identification of uncultured organisms but it can contribute to the noise observed in the mapping data, and thus to the need for stringent filtering.

As it relies on a learning step, the application of machine learning-based filters we have developed is limited to the training conditions. In our case, usage of random forest-based models to filter genomes should be limited to ecosystems with a similar complexity, sequenced with short reads at a depth included in the range used for the training and mapped to a catalog similar to representative genomes of UHGG in terms of completeness, intra and inter species diversity. More generally, readers should keep in mind that shallow shotgun metagenomics can be considered only for ecosystems for which (nearly) exhaustive reference databases exist, allows to assemble neither genes nor genomes, and can only produce coarse-grained functional analyses.

In addition to our analyses, we provide a pipeline that produces reliable taxonomic profiles at every sequencing depth, as an end-to-end solution for the analysis of shallow shotgun metagenomics data from human gut microbiota samples (see Availability of data and materials).

Conclusions

Our results show that (1) one needs to carefully filter taxonomic profiles retrieved from shallow shotgun metagenomics data to have trustworthy and interpretable results, (2) resulting taxonomic profiles are valid but limited to the most abundant taxa at low sequencing depths, i.e. taxa with relative abundance greater than 2.10 -4 at 1M reads/sample, and greater than 4.10 -4 at 500K reads/sample, and (3) shallow shotgun metagenomics allow a very good recovery of the structure of a data set, and constitutes a suitable approach to perform diagnosis-like classification of patients. It can be profitable in many gut microbiota-related clinical research projects to use shallow shotgun metagenomics. In comparison with metabarcoding, it allows to identify taxa down to species level and prevent biases related to targeted amplification. In comparison with deep shotgun sequencing, the information loss at cohort level is limited, thus it is certainly profitable to reduce the sequencing depth and favour number of samples to be analyzed to produce trustworthy biological knowledge.

Declarations

Figure 1 Simulations results: (A) Unambiguous fraction covered (uFC) and unambiguous read counts (uRC) of genomes present in the simulated profiles (TPs, blue points) or absent (FPs, red points). (B) Distribution of FN species according to their relative abundances in the simulated profiles, using RF-based filters with a 0.1 FDR on the testing set of cross validation. The dotted line represents the distribution of all simulated species. Tables Table 1 Classification of mapping hits in present and spuriously identified species: area under ROC curves and false negative rates when threshold is set to tolerate 0.1 FDR. For machine learning-based methods, these measures are split into training and testing sets, using a 4-fold cross validation.

Figure 2

 2 Figure 2 Comparison across sequencing depths for samples from the 3 data sets considered: richness observed, Shannon diversity and Bray-Curtis distance between subsampled data and reference (full depth data) using RF-based filtering (A, B and C respectively) and basic filtering (D, E and F respectively).

Figure 3

 3 Figure 3 Differences between patients groups in different studies: significance of inter-group difference regarding Shannon diversity index for α-diversity (A), PERMANOVA analysis for β-diversity (B). AUC corresponding to random forest classification (C) was performed in Loomba-2017 and Qin-2014, with a split between discovery and validation cohorts in Qin-2014 as performed on the original paper of this study. Error bars represent standard deviation for 10 repetitions of training process in the classification model.

Figures Figure 1

 1 Figures

Figure 2 Comparison

 2 Figure 2

Acknowledgements

We are grateful to the INRAE MIGALE bioinformatics facility (MIGALE, INRAE, 2020. Migale bioinformatics Facility, doi: 10.15454/1.5572390655343293E12) for providing help, computing and storage resources, as well as CRCM's DISC platform.

Funding

This work was financially supported by ANRT and Laboratoire Alphabio through a CIFRE scholarship (2019-0346).

Availability of data and materials

All data and scripts used in this article are available at https://forgemia.inra.fr/benoit.goutorbe/shallow-shotgun-metagenomics.git.

Abbreviations WGS : whole genome sequencing SSM : shallow shotgun metagenomics MAGs : metagenome-based assembled genomes UHGG : unified human gut genome [START_REF] Almeida | A unified catalog of 204,938 reference genomes from the human gut microbiome[END_REF] RC : read counts (mapped to a given genome) uRC : unambiguous read counts (mapped to a given genome) FC : fraction covered by all reads (for a given genome) uFC : fraction covered by unambiguous reads (for a given genome) RF : random forest

Ethics approval and consent to participate Not relevant

Competing interests

The authors declare that they have no competing interests.

Consent for publication

Not relevant

Authors' contributions BG, ALA, MM and SS conceived the ideas and designed the methodology; BG performed the simulations, analysed the data and developed the bioinformatics pipeline; BG led the writing of the manuscript with supervision from ALA, MM, AP, GB and SS. All the authors contributed critically to the drafts and gave final approval for publication.

Author details 1 Université Paris-Saclay, INRAE, MaIAGE, 78350, Jouy-en-Josas, France. 2 Clinical Research and R&D Department, Laboratoire Européen Alphabio, 13003 Marseille, France. 3 Cibi plateform, Centre de Recherche en Cancérologie de Marseille (CRCM), Institut Paoli-Calmettes, Aix-Marseille Université U105, Inserm U1068, CNRS UMR7258, 13009, Marseille, France.

Supplementary Files

This is a list of supplementary les associated with this preprint. Click to download. additional le1.xlsx listofadditional les.pdf


M M Larijani

A Navinrooz

F Le Normand

The bias-assisted HF CVD nucleation of diamond: Investigations on the substrate temperature and the filaments location

Keywords: Diamond, Bias and hot-filaments-enhanced chemical vapour deposition, Nucleation, Scanning electron microscopy, Raman spectroscopy

The bias-assisted HF CVD nucleation of diamond:

Investigations on the substrate temperature and the filaments location

Introduction

In the early times of the CVD diamond synthesis, one single step was carried out to grow a diamond film by the hot filaments CVD process (HFCVD), with previous exsitu treatments (by either carbon scratching, or seeding, or coating, or carburization, ...) to enhance nucleation [START_REF] Liu | Diamond Chemical Vapor Deposition: Nucleation and Early Growth Stages[END_REF]. Thus nucleation remained vanishingly little on mirrorpolished surfaces [START_REF] Constant | [END_REF]3]. However, it has been found, in a similar way to the microwave-assisted CVD process, that the application of a negative bias relative to the substrate before the growth in a short time strongly promotes the diamond nucleation, resulting in thin diamond films with even sometimes epitaxial orientation of the grains [4]. The deposition conditions such as substrate temperature, reactant composition, reactor pressure, etc., influence strongly the nucleation density, and the optimal conditions for the actual growth are not necessarily adequate for this initial nucleation step. Moreover, a balance must be found between the density and the quality of the diamond nuclei during the biasing step [5,6]. In a previous paper [7],w e described a new hot filament CVD set-up with direct current plasma assistance (DC HF CVD) for diamond nucleation. Activation of the gas phase constituents was obtained by a stable glow discharge between two grid electrodes coupled with two sets of parallel hot filaments located in-between and polarized at the corresponding plasma potential. In addition the sample is negatively biased with a small 10-15 V extraction potential with respect to the cathode grid. Such design allows creating a high density of both ions and radicals that are extracted and focused onto the surface of the sample. In this paper the influence of important nucleation parameters such as the distance between the filaments and the substrate and the substrate temperature is reported. While temperature is generally a key parameter of the nucleation by increasing the surface diffusion of atomic species, it was reported that the filaments-to-substrate distance might govern the energy and the concentration of the activated or ionized species [8]. Moreover, a direct comparison of the influence of these two parameters in the same experimental set-up might provide more insight into the true nucleation mechanism.

Experimental conditions

Si(100) samples, p-doped at 1-3 V cm and of size 8 Â 8 Â 0.5 mm 3 , were preliminary cleaned ultrasonically 5 min, consecutively into acetone and ethanol, then dipped into HF (4%) for 2 min. They were introduced into the CVD reactor (Fig. 1). The diamond deposition had been carried out in three consecutive steps: (i) surface etching by H + plasma, (ii) nucleation and (iii) growth. In the case of the etching and nucleation steps, a plasma was ignited both between the cathode and the anode grids (primary plasma: bias1 + bias2, Fig. 1) and between the cathode and the negatively biased substrate (extraction plasma: bias3, Fig. 1). The filaments were kept at the potential of the plasma (bias1, Fig. 1). The parameters related to the etching and growth steps were kept identical for all the samples. The substrate was independently heated by an infrared lamp located on the rear of the sample and controlled by a thermocouple. More details on the experimental set-up and the procedure were reported in Ref. [7]. Both to visualize images of the nuclei and to obtain a measurable Raman signal, the nucleation step was followed by a 1-h growth step. Table 1 summarizes these experimental conditions. The diamond films were characterized by Raman spectroscopy using the HeNe (632.8 nm) laser excitation source and a Renishaw set-up. A microscope objective was used to focus the laser onto a spot of approximately 1 Ami n diameter, and to collect the scattered light that then passed Fig. 1. Schematics of the DC HF CVD system. The distance d n between the substrate and the bottom couple of filaments is defined. The primary plasma between the cathode and the anode and the extraction plasma above the substrate are displayed in full. The parameters range of substrate temperature T and distance cathode -substrate are indicated in italics, while the most common used values are indicated in italic bold. V e and I e are the extraction voltage and current intensity, respectively (Fig. 1).

through the spectrometer onto a CCD detector. SEM investigations were performed with a JEOL JSM 840 working at 5 -10 kV.

Results

Nucleation temperature

The temperature of the nucleation was varied within the range 600 -850 -C. The Fig. 2 shows the Raman spectra and the SEM images at 600 and 850 -C, respectively. The SEM images illustrate that after 15 min of nucleation and 1 h of growth, the whole surface is covered by a diamond film. The film morphology is comparable for the whole range of nucleation temperature, with a b110À orientation. Although the grain size slightly increases with temperature, there are no appreciable differences in nucleation density. Raman spectroscopy is a convenient probe to evaluate the quality of the diamond films by recording the variations of the vibration frequency of the carbon -carbon bond. After subtraction of the luminescence background, the Raman spectra of Fig. 2 were decomposed into a lorentzian sharp peak D at 1332 cm À 1 characteristic of diamond and additional wide gaussian contributions. In all cases the Raman spectrum shows contributions centred at ¨1200, 1300, 1350, 1450 and 1560 cm À 1 [9]. The shoulders at ¨1200 and 1300 cm À 1 can be attributed to nanodiamond particles. Contributions at 1350 and 1560 cm À 1 are characteristic bands of sp 2 carbon clusters, the later one being due to the allowed mode E 2g and the former one is attributed to a breakdown in the selection rules due to size or defects-induced effects. Their relative contributions increase with nucleation temperature, as does the wide band centred around 1450 cm À 1 attributed to amorphous or DLC carbon with carbon in either threefold or fourfold coordination states.

Substrate-to-filaments distance

The substrate-to-filaments distance d n was varied from 4 to 8 mm, with the same set of filaments. Diamond particle size and density strongly change with the location of the filaments relative to the substrate. Thus the density decreases with increasing distance from a highly dense film at d n = 4 mm to a partial covering at d n = 8 mm, respectively, while the mean grain size increases. As shown in the Fig. 3 increase as the particle size decreases [10].A td n = 8 mm, the Raman intensity of diamond peak decreases significantly, in agreement with a significant decrease of nucleation density. Moreover a nanodiamond contribution around 1130 cm À 1 is quite noticeable.

Discussion

The nucleation density is quite sensitive to the location of the filaments but not to the temperature in the DC HF CVD process. By contrast film quality is little affected by both parameters, although a too little filaments-to-substrate distance leads to a decrease of the diamond quality and smaller size grain. Thus a distance around 5 -7 mm is believed to yield an optimized balance between a high nucleation density and a good quality of the diamond deposited. In a classical nucleation theory the nuclei rate formation is governed by the surface diffusion and formation of nuclei through surface collisions. Thus nucleation is strongly dependent on the surface temperature as surface diffusion is a thermally-activated process. By contrast in an ion-activated nucleation, enhanced diffusion is promoted by energetic ions impinging the surface. As the ion energy and concentration distributions are strongly dependent on the filaments location [8,11], it is expected that the nucleation formation rate will be now strongly modified with both the location and the thermal characteristics of the filaments. The filaments also play an important role to stabilize the primary plasma at a low polarisation voltage. Thus, with the presence of the filaments, the primary plasma is ignited around 200 V, much lower than without filaments [7]. This can be explained by the intense electron emission from the heated filaments that promotes the ionization reactions around the filaments. At low substrate-to-filaments distance we observe an asymmetric and broad diamond Raman peak with intense nanodiamond and amorphous carbon contributions. This is the indication of an intense secondary nucleation due to a high ion concentration above the substrate. At large filaments-tosubstrate distance reversely, the ion concentration decreases as a result of decreasing electron concentration and the nucleation rate decreases. Thus an optimization of the distance exists due to the spatially variable ion concentration in the extraction potential above the substrate.

The results provide evidence for an ion-enhanced nucleation rate. However two ion-enhanced nucleation processes were invoked in the literature [4]: Either surface collisions or bulk collisions through shallow implantation into the substrate, a process that is likely with energy of the ions larger some tenths of eV. The proportion of each of these processes depends on the kinetic energy of the ionic species impinging the surface. Thus it will also depend on the location of the filaments, all other parameters being otherwise constant. Further analyses are therefore required to discriminate between these two explanations and thus to elucidate the subtle behaviour of the filaments in the nucleation step of diamond.

Conclusion

The influence of substrate temperature and filaments-tosubstrate distance was studied during the diamond nuclea-tion step in a direct current and hot filaments-assisted CVD process. Whereas the nucleation temperature in the range of 600 -850 -C has negligible influence on the diamond nucleation density and film quality, the nucleation density is closely related to the distance between the filaments and the substrate. At low distance, secondary nucleation with nanodiamond formation and sp 2 -hybridized carbon are formed in addition to diamond. This is explained by an intense ion concentration above the surface. Reversely, at large distance, poor nucleation occurs due to decreasing ion concentration. It is thus expected an optimization of the distance to get a high nucleation density with a good film quality. However, not only the concentration, but also the energy of the ions, govern the diamond nucleation mechanism. Therefore combined surface studies and gas phase studies are required to elucidate the subtle role of the filaments in the nucleation step.

Fig. 2 .

 2 Fig. 2. Raman spectra and SEM images for different nucleation temperatures: (a) 600 -C, (b) 850 -C.

Fig. 3 .

 3 Fig. 3. Raman spectra and SEM images corresponding to substrate-filaments distance d n = 4 (I), 6 (II) and 8 mm (III).

Table 1

 1

	Experimental parameters during the diamond films deposition steps
	Step	T	P	Total flow rate	Distance cathode -	Filament temperature	Time	[CH 4 ]	V e	I e
		(-C)	(mbar)	(CH 4 +H 2 ) (sccm)	substrate (mm) d n	(-C)	(min)	(%)	(V)	(mA)
	Etching	700	15	100	4	1900	2	0	20	0.2
	Nucleation	600 -700 -850	15	100	2-4-6	1900	15	4	15	0.8
	Growth	850	30	200	6	2100	60	1	0	0


Joao Oliveira
email: joao.oliveira@irt-saintexupery.com

Bernardo Cougo
email: bernardo.cogo@irt-saintexupery.com

Fabio Coccetti
email: fabio.coccetti@irt-saintexupery.com

Stéphane Azzopardi
email: stephane.azzopardi@safrangroup.com

Hervé Morel
email: herve.morel@insa-lyon.fr

New Methodology for Defining Integration Limits Used for Switching Energy Computation in Power Devices

Keywords: Switching losses, SiC MOSFET, Hard switching, Calculation method, Wide bandgap devices

Switching losses of a power converter is a relevant factor that contributes to total losses. For building an optimized design, it is essential to accurately estimate switching energy in order to define certain parameters such as maximum operating frequency and cooling system. Power SiC MOSFETs have high switchingspeed capability, thus the voltage and current alignment should be performed carefully during classical dynamic characterization methods. In addition, the definition of integration limits of instantaneous power during switching events can generate erroneous results. In this paper, a new methodology for switching energy computation is proposed and validated for different manufacturers. This method

is less sensitive to oscillations on voltage and current waveforms with respect to other classical methods.

I. INTRODUCTION

Accurate switching energy measurement is essential to design power converters, especially those operating in high switching frequency. In addition, cooling system volume is directly impacted by the amount of dissipated energy [START_REF] Christen | Analytical switching loss modeling based on datasheet parameters for mosfet s in a half-bridge[END_REF]. Energy stored in the output capacitance of power SiC MOSFETs is much lower than that for MOSFETs based on Si material, which yields lower switching losses and higher efficiency in hard switching applications. However, this is also a source of spurious measurements of switching losses due to the resonance effect between capacitance and inductance elements associated to the test bench, and the power device [START_REF] Cougo | Accurate switching energy estimation of wide bandgap devices used in converters for aircraft applications[END_REF]. High frequency oscillations in the device, current and voltage, are found and become hard to separate the reactive power and the real dissipated power on the device under test. Classical methods for defining integration limits of device's instantaneous power have short integration duration, therefore they are often sensitive to these oscillations.

Different criteria act on limits integration definition of power signal. Reliability tests performed on SiC devices, especially on gate level, can generate voltage and current waveforms distortion during switching events [START_REF] Kozak | Degradation of sic mosfets under high-bias switching events[END_REF]- [START_REF] Witcher | Methodology for switching characterization of power devices and modules[END_REF]. During most of these tests, switching energy evolution is monitored according to an applied stress, thus an accurate measurement should be ensured. This can lead to erroneous switching energy measurement, depending on the considered initial and final switching instants.

This paper proposes a new methodology for defining limits integration used to switching energy computation. This method is less sensitive to oscillations on voltage and current waveforms with respect to other classical methods. Thus, compared to classical methods, the proposed methodology is more adapted to compute switching energy after aging tests. To validate this approach, switching energy is calculated before and after the reliability test GSS (Gate Switching Stress) [START_REF] Jedec | JEP195: Guideline for Evaluating Gate Switching Instability of Silicon Carbide Metal-Oxide-Semiconductor Devices for Power Electronic Conversion[END_REF].

II. SWITCHING ENERGY COMPUTATION

The DPT (Double Pulse Test) is widely used to analyze switching behavior of power devices, due to its simple circuit implementation [START_REF] Zhang | Methodology for wide band-gap device dynamic characterization[END_REF], [START_REF] Minasian | Power mosfet dynamic large-signal model[END_REF]. The main focus of the DPT is to characterize the dynamic performance, and performs switching energy measurements. The turn-ON switching energy (E on ) and turn-OFF switching energy (E of f ) are calculated by the integral of the drainsource voltage (V DS ) and drain-source current (I DS ) product, over on-transient (t on ) and off-transient t of f , respectively [START_REF] Jedec | JEP187: Guidelines for Representing Switching Losses of SIC MOSFETs in Datasheets[END_REF]. The key points for accurate E on and E of f calculation when DPT bench is used are: switching V-I timing alignment and integration limit definition of the power switching. Thereby, inaccurate measurements of switching energy can be associated to probe delays and bandwidth of instruments, and also to erroneous integration limit definition on switching power signal [START_REF] Zhang | Methodology for wide band-gap device dynamic characterization[END_REF]. WBG devices are more sensitive than Si-based devices due to their high-speed transients. Thus, even some nanoseconds of misalignment between voltage (V DS ) and current (I DS ) waveforms can generate significant error.

Delay compensations are widely discussed in the literature. In [START_REF] Chen | Characterization and modeling of high-switchingspeed behavior of sic active devices[END_REF], the deskew is performed by using a common square waveform generated by a function generator. One probe is used as baseline, and the other (voltage or current probe) is adjusted to align its rising or falling edge. DPT circuit is modified in [START_REF] Callanan | Sic mosfet double pulse fixture[END_REF] by removing the inductor and replacing the diode with a low inductance resistance. This allows for performing the deskew between voltage and current probes.

The total switching energy (E SW total = E on +E of f ) measured in DPT can be used to estimate switching losses. However, it is known that E on and E of f are not exactly the energy lost during turn-ON and turn-OFF, mainly because of the energy stored in the semiconductor device output capacitance (C OSS ). However, the integration over a switching period of the power corresponds well to the switching losses. This is due to the fact that over a cycle, the stored energy at the beginning and at the end of the cycle are the same, as the definition of a cycling system [START_REF] De Oliveira | A methodology for designing sic and gan device based converters for automotive applications[END_REF].

Based on a half-bridge power module, reference [START_REF] Zhang | Methodology for wide band-gap device dynamic characterization[END_REF] proposes calculating the difference between the input energy supplied by a DC capacitor and the output energy stored in a load inductor. However, this method takes into account a system without resistive losses, that can imply to an erroneous total switching energy.

A Modified Opposition Method was proposed in [START_REF] Cougo | Accurate switching energy estimation of wide bandgap devices used in converters for aircraft applications[END_REF] as a non-intrusive method for switching energy computation. It consists of an indirect measurement, in which four similar devices are tested in a full-bridge configuration through a test bench with accurately characterized elements. This method requires an accurate estimation of dissipated energy associated to circuit elements.

For this paper, an alternative to DPT bench is used. This new test bench allows to reduce total time the DC bus voltage is applied to the DUT, consequently it reduces the stress on DUT during switching characterization. More details will be presented in Section III.

A. Integration Limit Criteria

The definition of start and end points used to define integration limits on power signal (V DS I DS product) is a major challenge associated to DPT analysis [START_REF] New | On the estimation of switching loss in Wide Band-Gap power electronics systems[END_REF]. Depending on the considered initial and final switching instants, different energies can be calculated. In Table I, three methods from manufacturer application notes are presented for integration limit definition applied to power SiC MOSFETs.

10% final VGS 10% initial VDS 90% initial VGS 90% final VDS Wolfspeed [15] (Method B) 10% final IDS 10% initial VDS 10% final VDS 10% initial IDS Infineon [16] (Method C) 10% final VGS 2% initial VDS 90% initial VGS 2% initial IDS
Guidelines to select the integration limits cannot be universally applied for all SiC MOSFETs. A device with low rating current can have strong oscillations on gate input, therefore an integration window defined by use of gate-source voltage can generate spurious results, as it will be shown in Section V.

III. TEST BENCH FOR SWITCHING ENERGY MEASUREMENTS

The DST (Double Source Test) bench is an alternative to DPT that allows to stress a power device through high voltage/current switching. In Fig. 1 and Fig. 2, the electrical schematic and the used test bench are presented, respectively. With this test bench, it is possible to define independent applied values for the voltage and current on the DUT [START_REF] Garrab | On the extraction of pin diode design parameters for validation of integrated power converter design[END_REF]. In this way, the main goal is to perform switching characterization, i.e. switching energy measurement, with the possibility of varying several parameters, such as applied voltage/current (up to 1000 V and 100 A), switching speed (by changing gate resistance R G ), negative/positive applied gate voltages, etc. Signal generators are used to control the switches of the pulsed current source and the DUT. Switching signals are led to each gate driver boards by means of optical fiber. An auxiliary power supply is used to generate two different voltage levels: 12 V (for driver PCB) and 5 V (for optical isolator PCB). Two power supplies (not seen in Fig. 2) are used to define the required voltage (V A ) and current I L , which depends on (V F ). For the experimental tests shown in section V, three 1200 V SiC MOSFETs from different manufacturers were used. In Fig. 3, experimental waveforms of voltage and current on a SiC MOSFET can be seen. These signals were extracted from the DST bench and will be used as reference to present the new methodology for switching energy computation.

IV. NEW PROPOSED METHODOLOGY

Classical measurement precautions have to be taken in order to ensure large bandwidth by using suitable instruments for high speed signals. Moreover, the delay compensation between voltage and current probes is performed based on the propagation speed of traveling waves through the transmission line [START_REF] Oliveira | Switching behavior and comparison of wide bandgap devices for automotive applications[END_REF], [START_REF] Pérez | Electromagnétisme: fondements et applications avec 300 exercices et problèmes résolus[END_REF]. Integration intervals, as shown in Table I, are usually defined based on drain-source voltage and current waveforms. This can lead to some inconsistencies when high oscillations are observed due to the fast switching transients.

In addition, common-mode disturbances may cause abnormal effects on waveforms. Thus, the integration limit definition is very sensitive to all elements around the DUT.

The new proposed methodology, the Plateau Limits based Energy Computation (PLEC), is based on the switching power signal, as seen in Fig. 4. A plateau is defined by a flat part of a signal. It is described as a signal part included in a rectangle whose main parameters are a minimum width (duration) and a maximum height (fluctuation). The integration period of switching power signal is specified using plateau definition. The vertical fluctuation is set to 5% of the peak power in OFF or ON state. This parameter is experimentally validated by analyzing the primitive integral of power signal. It is important that the plateau duration is compatible with switching speed of the device, otherwise conduction losses will be taken into account in the switching energy computation.

Integration interval is defined by two points that correspond to two plateaus of 5% of the peak power. Inside the plateau, any integration interval produces a very similar energy. Thereby, this method is less sensitive to disturbances on voltage and current waveforms. Final integration point of the boundary should respect the plateau, where oscillations are nearly constant. The primitive integral of power signal corresponds to the switching energy. It was calculated for turn-OFF and turn-ON events, as can be seen in Fig. 5 and Fig. 6, respectively.

During turn-OFF transient, the oscillations are The number of positive areas (gray ellipses) is also a reference.

mainly due to a resonance between common-mode capacitance and output capacitance of DUT (C OSS ), and power loop inductance. This corresponds to an oscillating reactive energy stored essentially on the DUT, and it can cause erroneous measurements whether a short integration boundary is considered. After the plateau, it is possible to see a slight increase of the switching energy related to this phenomenon, as seen in Fig. 5.

For turn-ON transient, the parasitic capacitance discharging increases the power peak. In addition, at the end of the switching event, fewer oscillations are detected. After the plateau, the primitive integral has a negative slope due to a negative switching power, as seen Fig. 6. This is mainly related to an inaccurate measurement of the drop voltage on DUT during ON state, which can be compensated on signal post-processing.

PLEC method is based on power signal integral, and the final limit it is defined when minimum fluctuation of plateau is observed. This ensure that only switching losses is considered for computation. By definition of a plateau, a small change of a integration boundary inside a plateau do not change significantly the value of the computed switching energy.

V. COMPARISON OF THE DIFFERENT METHODS

PLEC methodology is compared to the methods mentioned in Table I. In order to have a better understanding of the difference on integration limits on power signal of a SiC MOSFET device, it is plotted, from experimental measurements, the drain-source voltage (V DS ) and current (I DS ), and gate-source voltage (V GS ), with the integration limits for three classical methods and the PLEC method. These results are shown in Fig. 7 SiC and Fig. 8. It is important to mention that the measured V GS is not the real signal at die level, but the voltage at the TO-247 terminals, which has higher oscillations than the real V GS at die level [START_REF] De Oliveira | A methodology for designing sic and gan device based converters for automotive applications[END_REF]. Thereby, methods for integration limits definition based on V GS signal can generate spurious results. Oscillations seen in Fig. 8 are mainly due to the strong coupling between power and gate loops, given that the used device has a 3-leads package, and also to common-mode current introduced by probe. The use of optically isolated probe could reduce these oscillations. The instantaneous power curve is generated after delay compensation performed between voltage and current probes. The integration limits on switching power according to each discussed method are depicted in Fig. 9, and Fig. 10. Comparison between the methods found in the literature, limited to 3 cases in this proposed study, and the PLEC method shows a large difference regarding the integration interval. A short integration period can be highly sensitive to oscillations due to the test bench or measurement instruments used during the tests. The insensitivity to oscillations is essential when it comes to aging tests that can modify harmonics on voltage and current waveforms. The device can have the same behavior, thereby presenting nearly the same total switching energy, but in a condition with sightly different oscillations.

Comparison between the mentioned methods, and the PLEC method for switching energy computation are presented in Table II, Table III, and Table IV. Three different power SiC MOSFETs, regarding rating current and manufacturer, are considered. Method A presents the highest deviation regarding the other methods. This is essentially due to the fact that integration limits definition is related to both V GS and V DS signals.

For device 1, the difference between total switching energy calculated by the PLEC method and others is very low. This is due to fact that this device presents lower oscillations among tested components. Device 2 presents the highest deviation which is about 38% regarding the PLEC method, when compared to the average value of total switching energy computed using methods A, B, and C. For device 3, a deviation of about 28% is observed.

VI. METHOD VALIDATION AFTER AGING TEST

Switching energy is often used as monitoring parameter during aging tests. Thereby, accurate measurements are essential to detect possible degradations. GSS reliability test is applied on several samples of a power SiC 11 and Fig. 12, the results are presented and compared to fresh devices (without stress).

After GSS test at 125 • C, method A presents higher losses when compared to after GSS test at 200 • C, which seems to be curious. Methods B and C do not reveal same trending after stress period at high temperatures. They present negative variation for the stress at 125 • C, and positive variation at 200 • C. Considering PLEC method, the energy variation after GSS test at 125 • C and 200 • C is about 3% and 10%, respectively. Higher oscillations found after stress period can modify integration limits, and thus generate erroneous measurements of switching energy. Proposed methodology presents more accuracy because the limits definition on power signal is more strictly defined, given that it is based on peak power and plateau definition.

VII. CONCLUSION

A new methodology (PLEC) for integration limit definition for switching energy computation in SiC power devices was proposed. This new method is based on the primitive integral of the switching power (V DS I DS product). A plateau is used as criterion for defining initial and final integration limits, by ensuring low dispersion of switching energy if the limits remains inside the plateau. The proposed method presents less sensitivity to oscillations when compared to classical methods found in the literature. This is essential, for instance, to correctly estimate switching energy when aging tests are applied to the DUT that can modify switching waveforms.

The PLEC method was defined and compared to classical methods using experimental waveforms from characterization of three different devices. After stress test applied on one device, switching energy was computed using these four different methods. PLEC method presented more coherent results regarding GSS stress in different temperatures.

Fig. 1 .

 1 Fig. 1. Electrical schematic of double source test.

Fig. 2 .

 2 Fig. 2. DST bench used for switching characterization.

Fig. 3 .

 3 Fig. 3. Waveforms of drain-source voltage and drain-source current in an inductive-load switching test. Green areas highlight the switching transient fronts.

Fig. 4 .

 4 Fig. 4. Instantaneous switching power of a SiC device. The criteria used to define the new integration intervals define a plateau. Parameters: τ of f = 400 ns, τon = 40 ns, and τp = 40 ns (adjustable).The number of positive areas (gray ellipses) is also a reference.

Fig. 6 .

 6 Fig. 6. Primitive integral of switching power during turn-ON event. Parameters: τon = 40 ns.

Fig. 7 .Fig. 8 .

 78 Fig.7. Waveforms of drain-source voltage and drain-source current during turn-OFF and turn-ON events. Red arrows correspond to the integration limits for the methods A and B, and the final integration limit for the method C.

Fig. 9 .

 9 Fig. 9. Switching power during turn-OFF event. The red arrows correspond to the integration limits for the methods A, B, and C. The blue arrows correspond to PLEC method. Parameters: τ of f = 400 ns, and τp = 40 ns (adjustable).

Fig. 10 .

 10 Fig. 10. Switching power during turn-ON event. The red arrows correspond to the integration limits for the methods A, B, and C. The blue arrows correspond to PLEC method. Parameters: τon = 40 ns, and τp = 40 ns (adjustable).

  presented in this paper. The goal is to perform a comparison between PLEC and other methods, before and after aging phase. The total time of stress test is 4000 hours, and the devices are stressed under two different temperatures that are 125 • C and 200 • C. In Fig.

		TABLE IV
	COMPARISON OF SWITCHING ENERGY COMPUTED ACCORDING
	TO DIFFERENT METHODS FOR A POWER SIC MOSFET
		(1200V;4.7A).
		DEVICE 3
	Operation Point (400V;20A)	Method A Method B Method C PLEC
	Turn-ON Energy (µJ)	19.4	19.0	19.7	19.6
	Turn-OFF Energy (µJ)	37.1	67.2	67.4	86.9
	Total Switching Energy (µJ)	56.5	86.2	87.1	106.5
	Total Switching
	Energy variation	0.53	0.81	0.82	1.00
	regarding PLEC
	MOSFET, and total switching energy is computed based
	on different methods

  Fig. 11. Total switching energy according to PLEC, A, B and C methods. GSS test is performed with devices under 125 • C.

					Total Switching Energy
		300
			VDS = 400V					Switching energy variation:
			IDS = 20A					PLEC:	ΔEsw= +3%
		250	@1200-33A device				A:	ΔEsw= +70%
								B:	ΔEsw= -10%
	E OFF (μJ) ON +E	200 150 100		Fresh devices		C: After Gate Switching ΔEsw= -14%
							Stress with T Tv,j v, j =125 °C =125 °C
		50
		0
			PLEC	A	B	C	PLEC	A	B	C
			Methodology for Defining Integration Limits
					Total Switching Energy
		300
			VDS = 400V					Switching energy variation:
			IDS = 20A					PLEC:	ΔEsw= +10%
		250	@1200-33A device				A:	ΔEsw= +6%
								B:	ΔEsw= +2%
	E ON +E OFF (μJ)	100 150 200		Fresh devices			C:	ΔEsw= +3%
		50					After Gate Switching
							Stress with T Tv,j v, j =200 °C =200 °C
		0
			PLEC	A	B	C	PLEC	A	B	C
			Methodology for Defining Integration Limits

Fig. 12. Total switching energy according to PLEC, A, B and C methods. GSS test is performed with devices under 200 • C.

This work was carried out in the framework of the IRT Saint Exupéry's project SiCRET (SiC Reliability Evaluation for Transport). We acknowledge the financial support from the SICRET industrial and academic members and the financial support from the French National Research Agency (ANR). SiCRET Industrial members: Alstom, Alter Technology, Emotors, Liebherr, Nucletudes, Safran, SuperGrid Institute, Thales, Vitesco The authors would also like to thank Sylvain Jouanolle from the French MoD (DGA) for the valuable support