We present a novel approach that combines machine learning based interactive image segmentation using supervoxels with a clustering method for the automated identification of similarly colored images in large image sets which enables a guided reuse of interactively trained classifiers. Our approach solves the problem of deteriorated segmentation and quantification accuracy when reusing trained classifiers which is due to significant color variability prevalent and often unavoidable in biological and medical images. This increase in efficiency improves the suitability of interactive segmentation for larger image sets, enabling efficient quantification or the rapid generation of training data for deep learning with minimal effort. The presented methods are applicable for almost any image type and represent a useful tool for image analysis tasks in general.

The provided free software TiQuant makes the presented methods easily and readily usable and can be downloaded at tiquant.hoehme.com.

Guided interactive image segmentation using machine learning and color-based image set clustering

Over the last decades, image processing and analysis has become one of the key technologies in systems biology and medicine. The quantification of anatomical structures and dynamic processes in living systems is essential for understanding the complex underlying mechanisms and allows, i.a., the construction of spatio-temporal models that illuminate the interplay between architecture and function. Recently, deep learning significantly improved the performance of traditional image analysis in cases where imaging techniques provide large amounts of data. However, if only few images are available or qualified annotations are expensive to produce, the applicability of deep learning is still limited.

Introduction

Advances in imaging technology have led to a diversification of microscopy techniques enabling examination of a broad spectrum of biological questions and have thus contributed to the establishment of imaging and image analysis as one of the main pillars of bioscience. Optical microscopy in combination with (immuno)histochemical and immunofluorescent staining allows for visualization of tissue and physiological entities. Therefore, these techniques are widely used [START_REF] Hammad | Protocols for staining of bile canalicular and sinusoidal networks of human, mouse and pig livers, three-dimensional reconstruction and quantification of tissue microarchitecture by image processing and analysis 2[END_REF] in medical diagnosis, e.g., using histological sections of tissue specimens or biopsies, as a means to distinguish physiological and pathological tissue architectures, and (2) in medical research to study tissue microarchitecture, e.g., using 3D volumetric confocal microscopy [START_REF] Hammad | Protocols for staining of bile canalicular and sinusoidal networks of human, mouse and pig livers, three-dimensional reconstruction and quantification of tissue microarchitecture by image processing and analysis 2[END_REF][2][START_REF] Hoehme | Prediction and validation of cell alignment along microvessels as order principle to restore tissue architecture in liver regeneration[END_REF], and tissue-scale or intracellular processes, e.g., by time resolved two-photon microscopy [START_REF] Vartak | Intravital multi-modal flux analysis reveals the transport mechanism of bile acids through hepatic microconduits[END_REF]. In many use-cases, computer-aided image analysis is preceded by detection or pixel-accurate segmentation of objects of interest. Detection asks for the location of an object in the image, e.g., using bounding boxes, while segmentation tries to capture the precise shape of the object. The most basic segmentation technique, manual pixel-accurate labeling, is time consuming, inconsistent and practically infeasible in many use cases, such as e.g., blood vessel segmentation in volumetric images. A multitude of (semi-)automatic image processing methods have been proposed, including e.g., intensity-thresholding and morphological operators [START_REF] Hammad | Protocols for staining of bile canalicular and sinusoidal networks of human, mouse and pig livers, three-dimensional reconstruction and quantification of tissue microarchitecture by image processing and analysis 2[END_REF][2][START_REF] Hoehme | Prediction and validation of cell alignment along microvessels as order principle to restore tissue architecture in liver regeneration[END_REF], region-based methods [5] or deformable models [START_REF] Mcinerney | Deformable models in medical image analysis: a survey[END_REF]. These methods are usually tailored towards a specific image setup (i.a., image dimensionality, magnification, staining) and/or object class (e.g., nuclei, blood vessels, necrotic tissue) and as such are typically limited in their applicability to differing use cases. Additionally, method parametrization to cope with image variability is often challenging and time consuming. More recently, shallow and deep machine learning methods were successfully applied to pixel-accurate image segmentation [START_REF] Lucchi | Supervoxel-based segmentation of mitochondria in em image stacks with learned shape features[END_REF][START_REF] Ren | Learning a classification model for segmentation[END_REF][9][START_REF] Ciresan | Deep Neural Networks Segment Neuronal Membranes in Electron Microscopy Images[END_REF]. By learning object appearance from training examples these methods allow for incorporation of subtle expert knowledge, are less restricted regarding image or object class specifics and minimize parametrization complexity.

Deep learning proved to be an extremely powerful approach to improving performance across all image processing tasks [9, [START_REF] Krizhevsky | ImageNet Classification with Deep Convolutional Neural Networks[END_REF][START_REF] He | Mask R-CNN[END_REF][START_REF] Long | Fully convolutional networks for semantic segmentation[END_REF]. This is primarily due to its ability to operate directly on the input image and implicitly learn features of increasing complexity. The price of finding a suitable feature space in addition to the decision surface itself is a high demand for training data and long training times. By now there is a large variety of specialized network architectures for segmentation of biomedical images and pre-trained models for standard targets such as nuclei or cells that can be applied directly if the images to be segmented are similar to images in the training set, or that can be tuned using transfer learning thereby reducing the amount of training data needed. However, the range of imaging systems and segmentation targets in the biomedical field is diverse so that often no suitable pre-trained models are available, image counts are in many cases very low and ground truth is usually not readily available and expensive to generate, which poses significant challenges for the applicability of deep learning methods. Furthermore, the development but also the training of deep learning models requires significant expert knowledge and is therefore difficult to achieve for untrained personnel. In contrast, traditional, so called shallow, machine learning approaches such as random forests or support-vector machines, rely on handcrafted features. While this limits versatility, as engineered features need to represent the learning problem at hand, finding a decision surface in a restricted feature space greatly simplifies the learning task. Thereby, relatively few training examples and short training times are required, which makes these shallow methods eminently suitable for interactive segmentation.

A number of interactive segmentation approaches and software packages that utilize machine learning have been published in recent years. A supervoxel-based approach for segmentation of mitochondria in volumetric electron microscopy images was developed by [START_REF] Lucchi | Supervoxel-based segmentation of mitochondria in em image stacks with learned shape features[END_REF], though the published executable software and code encompassed only the supervoxel algorithm SLIC.

SuRVos [START_REF] Luengo | SuRVoS: Super-Region Volume Segmentation workbench[END_REF] is a software for interactive segmentation of 3D images using a supervoxel-based hierarchy, which lays a focus on segmentation of noisy, low-contrast electron microscopy datasets. The extendable software package Microscopy Image Browser [START_REF] Belevich | Microscopy Image Browser: A Platform for Segmentation and Analysis of Multidimensional Datasets[END_REF] allows for processing of multidimensional datasets and features a selection of conventional processing algorithms, several region selection methods aiding manual segmentation, and methods for semi-automatic segmentation, including a supervoxel-based classification approach. FastER [START_REF] Hilsenbeck | fastER: a user-friendly tool for ultrafast and robust cell segmentation in large-scale microscopy[END_REF] is designed specifically for cell segmentation in grayscale images, using features that are very fast to compute. The Fiji plugin Trainable Weka Segmentation [START_REF] Arganda-Carreras | Trainable Weka Segmentation: a machine learning tool for microscopy pixel classification[END_REF] is a tool for pixel classification in 2D and 3D images using a broad range of features. Ilastik [START_REF] Berg | ilastik: interactive machine learning for (bio)image analysis[END_REF] provides a streamlined user interface with workflows for i.a., image segmentation, object classification and tracking as well as a sophisticated on-demand back end enabling processing of images with up to five dimensions that are larger than available RAM. Its segmentation workflow employs a random forest classifier for pixel classification from local features such as color, edge-ness and texture at different scales.

Most of these interactive segmentation tools allow for reuse of trained classifiers on new, unseen images. However, due to the variability of image appearance, even for sets of images that were acquired following a standard protocol, classifier reuse is usually a trial-and-error procedure and quickly becomes cumbersome for larger image sets. One of the main contributing factors to image quality variability are systematic discrepancies of coloration that can be attributed to minor deviations in the image acquisition process (e.g., sample preparation procedure, imaging settings, condition of the imaged subject). Since color information is used in various ways as a feature by all mentioned segmentation tools that target colored images, degraded prediction accuracy is to be expected for images with variable coloration.

We propose an approach that combines (i) interactive image segmentation from sparse annotations with (ii) the guided reuse of thereby trained classifiers on unseen images to enable efficient batch processing. (i) We formulate interactive pixel-accurate segmentation as a machine learning problem working on supervoxels, which are connected, homogeneously colored groups of voxels, using random forests or support vector machines as classifiers. Dimensionality reduction through use of supervoxels, precomputation of supervoxel features and a convenient graphical user interface enable rapid, intuitive refinement of training annotations by iterative correction of classification errors or uncertainties. (ii) We introduce a color-based image clustering method that enables automated grouping of image sets into subsets of similarly colored images.

A corresponding number of so-called prototype images is identified which serve as eligible candidates for interactive training of classifiers for within-subset reuse. These new methods, that are applicable to two-and three-dimensional images, have been fully incorporated into our freely available image processing software TiQuant. The first version of TiQuant, released in 2014, included various processing methods specifically for liver tissue segmentation of 3D confocal micrographs, as well as the corresponding analysis functionality [START_REF] Hammad | Protocols for staining of bile canalicular and sinusoidal networks of human, mouse and pig livers, three-dimensional reconstruction and quantification of tissue microarchitecture by image processing and analysis 2[END_REF]2].

We evaluate the interactive image segmentation method as well as the color-based image clustering strategy to guide classifier reuse using a previously published image set consisting of 22 brightfield micrographs of mouse liver tissue with corresponding manual nuclei annotations [START_REF] Hoehme | Prediction and validation of cell alignment along microvessels as order principle to restore tissue architecture in liver regeneration[END_REF]. We show that the interactive approach outperforms a human annotator and a comparable state-of-the-art interactive segmentation software and that limiting classifier reuse to similarly colored images significantly improves performance compared to reuse on images of differing coloration, yielding results close to the level of a human annotator.

Method

The general workflow of segmenting an image with our supervoxel-based approach is split into a preprocessing step and the interactive training, prediction and segmentation steps, that require direct user intervention and must be iteratively repeated for refinement as needed to achieve a segmentation of desired quality (see Fig. 1A).

In the preprocessing step, the image is partitioned into similarly sized supervoxels and descriptive features are computed for them. A single parameter, the target supervoxel size, must be adjusted by the user to ensure supervoxels fit the objects of interest. This step is computationally expensive compared to the interactive steps because it scales with the number of pixels and superpixels.

For example, on a modern computer, it takes ~1 minute to partition a 7.5k x 7.5k pixel image into 100k superpixels and six minutes to calculate all available features for them. However, this step is usually performed only once per image.

Subsequently, the user annotates exemplary fore-and background regions in the original image to generate a training database. A classifier is then fitted to the training data, and probability estimates for class membership are predicted for all supervoxels in the image. In the final step, a segmentation is generated based on these probability estimates, which can be refined through post-processing.

If segmentation quality is insufficient, the user can provide further training data, especially for regions that were poorly segmented. Thereby, high-quality segmentations can be generated quickly by iterative annotation refinement.

Trained classifiers can be used for prediction on unseen images, eliminating the need for producing training data for every image. To guide the process of identifying training images and matching candidate images for classifier reuse, the images to be processed are grouped into similarly colored subsets using a two-stage clustering approach (see Fig. 1B).

Supervoxel initialization

The term supervoxel (respectively superpixel in 2D) was introduced by Ren and Malik [START_REF] Ren | Learning a classification model for segmentation[END_REF], who proposed oversegmentation, i.e., the process of partitioning an image into contiguous regions so that the objects of interest are themselves subdivided into distinct regions, as a preprocessing step to reduce image complexity while preserving most of the structure necessary for segmentation at the scale of interest. Supervoxels group voxels into visually meaningful building blocks of more or less similar size and compactness, depending on the chosen generative approach and parametrization. They reduce data dimensionality with minimal information loss, thereby greatly reducing computational cost of subsequent processing steps, and they allow the computation of local features such as color histograms and texture descriptors. To date, a multitude of different supervoxel algorithms exist which can be categorized by their high-level approach, into e.g., graph-based, density-based and clustering-based algorithms [19].

We use the SLIC0 variant of the clustering-based algorithm simple linear iterative clustering (SLIC) [START_REF] Achanta | SLIC superpixels compared to state-of-the-art superpixel methods[END_REF], due to its comparatively strong performance regarding Boundary Recall and Undersegmentation Error [19], and its memory and runtime efficiency [START_REF] Achanta | SLIC superpixels compared to state-of-the-art superpixel methods[END_REF]. Its iterative nature enables straightforward runtime restriction and provides direct control over the number of generated superpixels. SLIC is an adaptation of k-means clustering with two main distinctions: i)

The search space is reduced to a region proportional to the superpixel size, yielding a complexity linear in the number of pixels and independent of superpixel count. ii) The weighted distance measure used combines color and spatial proximity, providing control over size and compactness of the resulting superpixels. The SLIC0 variant adaptively choses the superpixel compactness, thereby reducing the free parameters to be adjusted to the number of superpixels, or the superpixel size, respectively. The SLIC0 algorithm was implemented as an ITK (Insight Toolkit) filter and extended to work in three dimensions. Examples of superpixel oversegmentations can be seen in Fig. 2.

In brightfield and confocal microscopy, tissues or physiological entities are stained; accordingly, color, color gradients or boundaries, and local color patterns are the primary source of information for image interpretation. Therefore, the implemented supervoxel features comprise local and neighborhood color histograms, edge-ness over several spatial scales, as well as texture descriptors (detailed description in SI Appendix), yielding a feature vector with a maximum of 138 entries. Feature categories can be disabled by the user to speed up processing.

Interactive training data generation

For training data generation, TiQuant provides a graphical user interface that allows visualization of and interaction with images. The user can draw directly on the image in order to denote exemplary regions for the classes to be segmented. Supervoxels in the annotated regions are collected and their feature vector together with the annotated class are written into a training database. For examples of how training data annotation is presented in the software, see Fig. 2.

Supervoxel-Classification

Random Forests (RF) and Support Vector Machines (SVM) are supported as classifiers to learn the class membership of supervoxels given their feature vectors. We included both classifiers to choose from as they were found to be the two best-performing classifier families when evaluated on a diverse set of learning problems [START_REF] Fernández-Delgado | Do We Need Hundreds of Classifiers to Solve Real World Classification Problems?[END_REF]. Furthermore, both classifier families are able to handle non-linear data and have few parameters, however, RFs are faster in regard to training and execution and are therefore set as default.

The user-supplied training data, which is summarized in the training database, is split into a training set and a test set in a 70:30 ratio to allow evaluation of classifier performance (see Fig. 3).

Since training data can be relatively sparse, we use this skewed ratio to ensure a sufficiently large database for training, potentially at the expense of accuracy in assessing classifier performance. Prior to classifier training, hyperparameter optimization of the chosen RF or SVM classifier can optionally be performed to tune parameters that are not directly learned to the learning problem.

In order to limit execution time while retaining the explorative quality of an exhaustive Grid Search, the optimization is done using Random Search, which tests a fixed number of parameter settings sampled from given distributions [START_REF] Bergstra | Random Search for Hyper-Parameter Optimization[END_REF]. The search is performed on the training split with 5-fold cross-validation employing the balanced accuracy score [START_REF] Brodersen | The Balanced Accuracy and Its Posterior Distribution[END_REF] to evaluate classifier performance.

As an intermediate step, the (optionally) optimized classifier is trained on the training split and evaluated on the withheld test split to assess its performance. This allows the user to evaluate the suitability of the selected features and classifier, to compare different configurations and to determine whether additional training data is expected to improve performance.

Finally, to make efficient use of the sparse training data, the classifier is trained and calibrated on the entire user-supplied data corpus using the previously determined optimized parametrization.

The resulting trained classifier is then used to predict the class membership probability estimates of all supervoxels of the image. Exemplary probability maps are shown in Fig. 2. Classifier calibration is performed after training using 5-fold cross-validation to provide calibrated probabilistic estimates of class membership as output. In case of SVMs, Platt Scaling [START_REF] Platt | Probabilistic outputs for support vector machines and comparisons to regularized likelihood methods[END_REF] is used for this purpose, fitting an additional sigmoid function to map SVM scores to probabilities. RFs, on the other hand, provide probabilistic estimates per default, which can be optionally calibrated using Platt Scaling or the non-parametric Isotonic Regression approach [START_REF] Barlow | Statistical inference under order restrictions; the theory and application of isotonic regression 26[END_REF]. Empirical results show that SVMs and RFs are among the models that predict the best probabilities after calibration [26]. The Brier score is used to evaluate classifier performance during calibration, since it is a proper scoring rule that measures the accuracy of probabilistic predictions [START_REF] Gneiting | Strictly Proper Scoring Rules, Prediction, and Estimation[END_REF].

Segmentation

The probabilistic estimates are mapped to a binary class-membership by applying a threshold, which defaults to the value of 50%. Subsequently, the resulting segmentation can be postprocessed using three optional filters. The first post-processing step allows for removal of isolated foreground objects smaller than a specified minimal foreground object size, as well as filling of holes in foreground objects smaller than a minimal hole size. Additionally, a morphological closing operator followed by an opening operator can be applied to smooth surfaces. Specialized label set operators are used that work on the label image representation in which pixel values represent connected objects. These operators preserve the initial label assignment by not merging disconnected foreground objects [START_REF] Beare | Parallel Algorithms via Scaled Paraboloid Structuring Functions for Spatially-Variant and Label-Set Dilations and Erosions[END_REF]. Finally, the watershed algorithm may be applied in order to split wrongly connected objects such as nuclei [START_REF] Malpica | Applying watershed algorithms to the segmentation of clustered nuclei[END_REF]. Two exemplary segmentations are shown in Fig. 2.

Color-based image clustering

A common but severe problem hampering efficient batch processing of a large number of images is visual heterogeneity of images due to varying imaging settings and conditions, or systematic visual changes during the imaged process. One approach to solve this problem is image normalization that aims to visually equalize images prior to segmentation and analysis. There is a large number of normalization methods available today, typically aiming at normalizing the intensity distribution of images and/or improving the signal-to-noise ratio, such as variations of adaptive histogram equalization [START_REF] Pizer | Adaptive histogram equalization and its variations[END_REF], color deconvolution [START_REF] Ruifrok | Quantification of histochemical staining by color deconvolution[END_REF], or more recent deep learning based methods [START_REF] Litjens | A survey on deep learning in medical image analysis[END_REF]. However, the selection of a suitable method as well as its parameterization is non-trivial. Moreover, normalization methods can easily introduce systematic errors in subsequent image analysis if differences in visual appearance originally result from changes in tissue conditions to be analyzed. In this case, normalization could remove properties that actually should be quantified.

We propose a different method to overcome the problem of visual heterogeneity in microscopic images by automatically and robustly grouping images of a given dataset into visually homogeneous subsets. Thereby, we can not only guide the selection of training images, but also facilitate classifier reuse leading to minimized manual annotation effort. The underlying assumption is that minor procedural discrepancies in the image acquisition process introduce systematic changes of coloring. Possible causes comprise slight variations between sample preparation sessions (e.g., affecting staining penetration depth, thus color saturation), minor deviations in imaging settings between imaging sessions (e.g., affecting brightness and contrast), as well as differing conditions of the imaged subjects (e.g., healthy vs. impaired tissue).

In a first step images are analyzed for their dominant colors (see Fig. 3). This process, also known as palette design, is one of two phases of color quantization, an operation used for e.g., image compression. It has been shown that k-means clustering is an effective method for this task [START_REF] Kasuga | Color quantization using the fast K-means algorithm[END_REF][START_REF] Celebi | Effective initialization of k-means for color quantization[END_REF]. In order to identify a set of the kc most dominant colors of an image in RGB color space, each pixel's RGB color vector is interpreted as a data point in 3D space. Instead of starting with fully random cluster centers the k-means++ initialization scheme [START_REF] Arthur | k-means++: The advantages of careful seeding 36[END_REF] is used, as it has been demonstrated to improve effectiveness for this task [START_REF] Celebi | Effective initialization of k-means for color quantization[END_REF]. Elkan's algorithm [36], which is a faster variant of Lloyd's algorithm [START_REF] Lloyd | Least squares quantization in PCM[END_REF], using the triangular inequality to avoid many distance calculations, is employed to cluster the RGB color data points.

Next, the entries of the dominant color sets are sorted, since their initial order is based on color prevalence, i.e., the number of pixels assigned to a respective color cluster. This ordering, however, is susceptible to changes in image composition, so that, e.g., the size or number of physiological entities influences the rank of the color cluster(s) they are assigned to. Sorting, which is done component-wise on the RGB vector entries, allows comparing the dominant color sets of different images.

Finally, the image set is partitioned into subsets with similar dominant color sets. Our approach extends previous work, in which color moments [START_REF] Maheshwary | Retrieving Similar Image Using Color Moment Feature Detector and K-Means Clustering of Remote Sensing Images[END_REF] or histograms [START_REF] Malakar | Image clustering using color moments, histogram, edge and K-means clustering[END_REF] were used as image descriptors for k-means based image clustering. Each image is represented by its 3*kc dimensional sorted dominant color set. K-means clustering is applied to this set of data points yielding ki clusters, minimizing the within-cluster variances of the sorted dominant color sets. The images with the smallest Euclidean distance to their respective cluster center are recommended training images for their image cluster.

Experiments Validation

We validated our interactive segmentation and image set clustering method exemplarily on a set of 22 brightfield micrographs of mouse liver paraffin slices. These images were hand-labeled by a trained expert to automatically extract architectural and process parameters, which were used to parameterize a spatio-temporal tissue model of tissue regeneration after CCl4 intoxication [START_REF] Hoehme | Prediction and validation of cell alignment along microvessels as order principle to restore tissue architecture in liver regeneration[END_REF]. Tissue slices were immunostained for BrdU positive nuclei to visualize proliferation and typically show a single liver lobule centered by a central vein. Image data was obtained for a control (t=0) and at 7 time points after administration of CCl4, which causes necrotic tissue damage in the area around the central vein. Over the observed time period, this damage regenerates and the lost cells are replenished by cell proliferation of the surviving liver cells.

We reuse this image set to validate our method as it is representative for many image segmentation tasks in a biomedical context, which are characterized by a small total number of available images and heterogeneous image appearance caused by variations in sample preparation, imaging settings, or tissue conditions, or by various other possible technical problems like entrapped air, blurring, or imperfect staining application.

A re-examination of the original nuclei annotations revealed several inaccuracies and inconsistencies (Fig. 4A). Therefore, two trained experts manually revised all 22 images thoroughly to define a gold standard nuclei annotation. The annotations represent each nucleus as a 2D pixel coordinate. As a measure for segmentation accuracy, we use the F1 score, thereby considering both precision and recall. The number of true positives (tp), false positives (fp) and false negatives (fn) underlying the score are quantified based on an object-wise mapping of the gold standard to the respective segmentation. Specifically, a segmented nucleus is considered a tp if a gold standard nucleus (i.e., its 2D coordinate) lies within the segmented region. Accordingly, a segmented nucleus within which there is no gold standard nucleus is counted as a fp. A fn is registered if there is no segmented nucleus in the immediate neighbourhood (kernel of size 1) of a gold standard nucleus. Based on the F1 scores, we validate and compare the following methods: i) original manual annotation (Fig. 4A), ii) our superpixel-based interactive processing approach (Fig. 4B), iii) intra-and iv) inter-cluster reuse of interactively pre-trained classifiers as well as v) the state-of-the-art tool ilastik [START_REF] Berg | ilastik: interactive machine learning for (bio)image analysis[END_REF] as a reference point.

For validation of the superpixel-based approach, each image was partitioned into approximately 50k superpixels of target size 64 pixels (initially 8x8 pixels), and each superpixel was analyzed for all available features, comprising gradient magnitude, Laplacian of Gaussian, local and neighborhood color histograms as well as texture features.

In order to evaluate the suitability of our tool for interactive segmentation, appropriate training data was produced for each of the 22 images by annotation (following workflow Fig. 1A). Per image, a RF classifier was trained and calibrated after hyperparameter optimization was performed. The trained classifiers were applied to the respective image to predict class membership of all contained superpixels. The resulting segmentations were post-processed by removing objects smaller than three superpixels, smoothing of boundaries of segmented nuclei, and finally by applying the watershed algorithm to split up agglomerations of nuclei into individual objects.

Subsequently, to evaluate how well the trained classifiers generalize to unseen images and whether a restriction to images with similar coloration improves performance, we grouped the image set into ki=6 subsets. The number of clusters was determined experimentally using the Elbow method, a heuristic for determining the optimal number of clusters representing the majority of the variability in a dataset. The number of dominant color sets kc=4 was chosen manually to correspond with the different image compartments of relevant size (i.e., background, nuclei, and two tissue compartments for healthy tissue/necrotic lesion, or glutamine synthetase positive/negative zone, respectively). Per image subset the image with the smallest Euclidean distance to the image cluster center was selected as prototype image. The previously trained classifiers of these prototype images were applied (i) intra-cluster to all other images in the respective cluster and (ii) inter-cluster to all images not belonging to the respective cluster.

Segmentation post-processing was done as described before.

We compare our approach with the established image processing software ilastik. For the results to be comparable, we used ilastik's 'Pixel Classification + Object Classification' workflow. All 37 predefined features were selected for training. The training annotations created with our software were imported and an object size filter was used for segmentation post-processing.

Runtime evaluation

A detailed discussion of aspects of computational performance is given in the SI Appendix. The interactive processing steps training and prediction are independent of pixel count and linear in the number of supervoxels. For the segmentation step, runtime scales linear with image size and number of supervoxels. Exemplarily, execution times for an image with 100M pixels and 100k supervoxels are ~2min and 9min for supervoxel generation without and with feature generation.

Training database generation, prediction and segmentation were completed in 5sec, 12sec and 51sec, respectively. The analysis was done on an Intel Core i9-7900X with 10 logical cores at a 3.3 GHz clock speed with 64 GB RAM.

Discussion

Currently, interactive image segmentation often is the optimal image processing approach for applications where a detailed analysis of objects of interest is the aim, but only few images are available, and/or training annotations are not provided. However, these seemingly specific conditions are fulfilled in many cases of biological and medical image processing where often pixel-accurate segmentations of physiological entities are required based on two-or threedimensional images. This includes the counting and measuring of (sub-)cellular structures and vascular networks in tissue. Custom-made image processing solutions that target specific physiological entities imaged utilizing a specific technical setup are still a very common and widespread approach. Unfortunately, these solutions usually cause considerable development overhead compared to interactive solutions based on machine learning that can be easily adapted to new data through training. Moreover, custom-made solutions require profound expertise in image analysis techniques including knowledge of particular algorithms and the impact of parameter choices. In contrast, interactive image segmentation using machine learning from manual image annotations has a much flatter learning curve and greatly simplifies the process of quantifying the content of complex images enabling users to directly apply their expert knowledge of, for example, biological or medical entities in the analyzed images. The method described in this paper reduces the technical parameters virtually to zero. Only the size of the superpixels has to be set, however, it can be inferred from the size of the structures to be quantified.

Furthermore, interactive segmentation is also applicable for larger datasets consisting of many images, when no ground truth or appropriate pre-trained deep learning models are available. In this case, interactive image segmentation can be useful in producing training data for deep learning with significantly reduced manual effort and thus cost.

The presented supervoxel approach using the SLIC0 algorithm is especially suited as a basis for responsive, interactive segmentation not only because it is easy to use, but also technically due to its dimensionality reduction characteristic, minimizing processing times and memory consumption. The latter effect is particularly pronounced if the objects of interest are highly resolved. On the other hand, the advantage diminishes if object size approaches those of single pixels. For objects with fluent boundaries or low-resolution images, the alignment of supervoxels to object boundaries can be suboptimal and traditional pixel-level segmentation methods might provide higher segmentation accuracy especially in boundary areas.

The features implemented in the presented work are designed for application to color images and might be less effective in representing learning problems formulated based on grayscale or binary images, for example electron microscopy or x-ray, for which optimized features could be integrated.

Learning from annotations on single images and reusing these interactively trained classifiers for unseen images would further decrease manual effort. However, the downside of learning from such sparse annotations is a decreased generalizability which is especially deteriorating results in case of high color variations between the images to be quantified. We have demonstrated that limiting classifier reuse to similarly colored images significantly increases the average performance on a validation dataset, and thus improves the suitability of interactive image segmentation for processing moderately sized image sets.

Conclusion

We present a novel approach that combines machine learning based interactive image segmentation using supervoxels with a clustering method for the automated identification of similarly colored images in large image sets which enables a guided reuse of interactively trained classifiers. Our approach solves the problem of deteriorated segmentation and quantification accuracy when reusing trained classifiers which is due to significant color variability prevalent and often unavoidable in biological and medical images. We demonstrate that our interactive image segmentation approach achieves results superior to both manual annotations done by a human expert and also an exemplary popular interactive segmentation tool when dedicated training data is provided for each image (see Tab.1). We show that by limiting classifier reuse to automatically identified images of similar color properties, segmentation accuracy is significantly improved compared to a traditional, non-discriminative reuse. In addition, our approach greatly reduces the necessary training effort. This increase in efficiency facilitates the quantification of much larger numbers of images thereby improving interactive image analysis for recent technological advances like high throughput microscopic imaging or high-resolution video microscopy [START_REF] Ljosa | Annotated high-throughput microscopy image sets for validation[END_REF].

Additionally, the presented method can readily be used for a fast generation of training data for deep learning approaches. In summary, our approach opens up new fields of application for interactive image analysis especially in but not limited to a biological and medical context. The provided free software makes the presented methods easily and readily usable.

Fig. 1 :

 1 Fig.1: (A) Workflow of the supervoxel-based interactive image segmentation approach from a user perspective. The classifier training procedure is illustrated in detail in Fig. 3A. (B) Workflow of guided reuse of interactively trained classifiers: The set of n images is grouped into m subsets of similarly colored images using our color-based image clustering method. For each subset a prototype image is identified, which is then interactively segmented. The resulting m trained classifiers are reused on the remaining images in their respective subsets. The color-based clustering algorithm is depicted in detail in Fig. 3B.

Fig. 2 :

 2 Fig.2: Illustration of supervoxel-based image segmentation procedure on an image of a cheetah (left), and a 3D confocal micrograph of mouse liver tissue in which blood vessels were segmented (right). The capillaries in this imaging setup appear yellow-green and usually comprise an unstained lumen as well as few small, elongated nuclei colored blue, which typically belong to endothelial (sinusoidal) cells forming the capillary wall, immune cells residing in the lumen or other non-parenchymal cells. (1) Supervoxel outlines are shown in gray. (2) Training data for the background class is colored red in both instances, while the foreground class is colored blue in the left and white in the right example. Annotations were selected to represent the characteristic visual features of both classes. For example, annotations of capillaries comprise yellow-colored walls, lumen and enclosed cell nuclei, whereas the background annotations encompass unstained cytoplasm and nuclei of the parenchymal liver cells (blue), bile canaliculi (green) as well as regions close to capillaries to promote learning of exact boundaries. (3) Class membership probabilities in the prediction images are illustrated using a color mapping ranging from red (low probability of being foreground) over yellow to blue (high probability of being foreground). (4) The segmentation is visualized by a blue overlay on the left, and a yellow overlay on the right.

Furthermore, the training

  data is expected to be imbalanced, i.e., that classes are unevenly distributed. To account for this disparity, the initial training-test split, as well as splits performed during cross-validation, are done in a stratified fashion to preserve relative class frequencies. Moreover, in the training phase, the training samples are weighted, where the weight is inversely proportional to the class frequency, and appropriate scoring functions that are insensitive to imbalanced training data are used for classifier evaluation. The feature vectors computed during supervoxel initialization, described in detail in the SI, are normalized by subtracting the mean and scaling to unit variance independently for each feature component to ensure comparable feature scales.

Fig. 3 :

 3 Fig.3: (A) Flowchart of the classifier training procedure. The set of annotated supervoxels is split in a 70:30 ratio into a training and test split. The RF or SVM classifier is trained on the training split after optional hyperparameter optimization and evaluated on the test split to allow the user to assess and compare classifier performance. The classifier employed for final prediction uses the optimized parametrization and is trained and calibrated on the entire original set of annotated superpixels. Fig 1A illustrates the integration of the training procedure into the interactive segmentation workflow. (B) Flowchart of the color-based image clustering algorithm. For each image, a characteristic set of the kc most dominant RGB colors is determined by k-means clustering of the RGB pixel values using the k-means++ initialization scheme. To ensure comparability of these representative sets, their RGB entries are sorted by component-wise comparison of R/G/B values. Finally, k-means clustering is applied in the 3*kc-dimensional space of sorted dominant color sets to group the corresponding n images into ki subsets. For each of these subsets, a prototype image is identified, which is recommended for interactive segmentation. Fig 1B illustrates the integration of the colorbased image clustering algorithm into the workflow of guided reuse of interactively trained classifiers.

  The summarized results of the method validation are shown in Tab.1, analysis is shown in the SI Appendix Fig.S3. The best result was achieved with the fully interactive approach where dedicated training data was provided for each image.The original manual annotation produced comparable results to intracluster reuse of classifiers trained on prototype images, indicating that the restriction of classifier reuse to similarly colored images is able to produce human-grade results. However, as the increased variance for classifier reuse indicates, segmentation quality is less reliable. Reusing classifiers on images with differing coloration generally degrades segmentation quality greatly, as expected given the importance of color information. Fully interactive segmentation with ilastik, using dedicated training data for each image, performed slightly better than the pre-trained classifiers on intra-cluster images.

Fig. 4 :

 4 Fig.4: Visualization of error classes from comparison with gold standard (GS): (A) Comparison of GS with the original manual segmentation created by a single human annotator: Correctly segmented nuclei (true positives: bright green); Incorrectly segmented nuclei (false positives: yellow); Nuclei that were falsely not segmented (false negatives: red). (B) Comparison of GS with a segmentation generated by our interactive segmentation tool: Coloring as in (A). Nuclei outlines in (A) were generated from 2D pixel coordinates provided by the annotator, while outlines in (B) were generated directly from the nuclei segmentation.

Institute of Computer Science, Leipzig University, Haertelstr. 16-18, 04107, Leipzig, Germany;

IfADo -LeibnizResearch Centre for Working Environment and Human Factors, Ardeystrasse 67, Dortmund, Germany;

INRIA Paris & Sorbonne Université LJLL, 2 Rue Simone IFF, 75012, Paris, France; *Corresponding author


Muriel Fartoukh

Guillaume Voiriot

Laurent Guérin

Jean Damien Ricard

Alain Combes

Morgane Faure

Sarah Benghanem

Etienne De Montmollin

Yacine Tandjaoui-Lambiotte

Antoine Vieillard-Baron

Jean Damien Ricard

Eric Maury

Jean-Luc Diehl

Keyvan Razazi

Virginie Lemiale

Pierre Trouiller

Benjamin Planquette

Laurent Savale

Nicholas Heming

Jonathan Marey

Ap-Hp Clinical

Data Warehouse

Fabrice Carrat

Nathanael Lapidus

The Epicuflu_Aphp Group

Seasonal Burden of Severe Influenza Virus Infection in the Critically Ill Patients, Using the Assistance Publique-Hôpitaux de Paris Clinical Data Warehouse: A Pilot Study

Keywords: Epidemic, Influenza, Assistance Publique-Hôpitaux de Paris (AP-HP) clinical data warehouse, Critical care, Prognosis

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

Background

Each year, the characteristics of the flu epidemic are likely to evolve [START_REF] Petrova | The evolution of seasonal influenza viruses[END_REF][START_REF] Nelson | The evolution of epidemic influenza[END_REF][START_REF] Harper | Seasonal influenza in adults and children-diagnosis, treatment, chemoprophylaxis, and institutional outbreak management: clinical practice guidelines of the Infectious Diseases Society of America[END_REF][START_REF] Delabre | Immunity against influenza A(H1N1) infections is determined by age at the time of initial strain circulation[END_REF][START_REF] Lemaitre | Mortality burden of the 2009 A/H1N1 influenza pandemic in France: comparison to seasonal influenza and the A/H3N2 pandemic[END_REF], and may require specific recommendations from the health authorities. In France, the Open Access *Correspondence: muriel.fartoukh@aphp.fr 1 Service de Médecine Intensive Réanimation, Hôpital Tenon, Assistance Publique-Hôpitaux de Paris, Sorbonne Université, 75020 Paris, France Full list of author information is available at the end of the article coordination of the epidemiological and virological surveillance has been gradually structured under the auspices of the National Institute for Public Health (Santé Publique France, SPF) (http:// invs. sante publi quefr ance. fr). At the critical care level, the surveillance system has been developed since the 2009 influenza pandemic, with heterogeneous regional modalities of implementation. In Paris area, individual case report forms of probable or confirmed severe influenza illness are completed by a regional network of 17 sentinel adult medical intensive care units (ICU) and affiliated medical intermediate care wards (ICW) on a voluntary basis, and are sent to SPF every time a new case is diagnosed. Based on those reports, a regional feedback is weekly available, to provide detection and situational awareness regarding the most severe cases in the Paris area.

The relevance of a real-time computerized tool for monitoring and reporting severe cases and their impact on critical care services has been assessed during the H1N1 2009 pandemic in a study involving the REVA research network in connection with the French Intensive Care Society [START_REF] Mercat | Severe H1N1 2009 influenza infection in adults: the French experience[END_REF], by demonstrating the impact of severe cases on the workload and organization of ICUs. Other efforts have been conducted to improve the assessment of the influenza epidemic severity and its impact on critical care [START_REF] Boddington | Developing a system to estimate the severity of influenza infection in England: findings from a hospital-based surveillance system between 2010/2011 and 2014/2015[END_REF][START_REF] Ercole | Modelling the impact of an influenza A/H1N1 pandemic on critical care demand from early pathogenicity data: the case for sentinel reporting[END_REF]. As part of the implementation of a common clinical information system for all the centers of the Assistance Publique-Hôpitaux de Paris (AP-HP) hospital group, the healthcare data warehouse "Entrepôt de Données de Santé" (EDS AP-HP https:// eds. aphp. fr) was set up in 2015 to support non-interventional research and hospital management studies based on data collected during patients' stays at AP-HP. The EPIcuFLU_APHP research is a pilot multi-site epidemiological study of admissions in adult critical care units (medical ICU/ICW and respiratory ICW) for Influenza virus infection. The main objective was to assess the burden of the epidemic on critical care units, by describing the severity and outcomes of adult patients admitted to the ICUs/ICWs of the APHP network during the influenza season. The primary and secondary endpoints were the in-hospital mortality within 28 days of ICU/ICW admission with a diagnosis of influenza infection; ICU/ICW and hospital lengths of stay and in-hospital mortality rates, and the early prognostic factors associated with 28-day mortality, based on data available during the first 24 h of ICU/ICW admission.

Materials and methods

Study design and population

The research was conducted during the 2017-2018 influenza epidemic in France [START_REF]Equipes de surveillance de la grippe B-S Sibylle. Influenza activity in France, season 2017-2018[END_REF], from November 1st 2017 to May 31 2018, in the medical adult ICUs/ICWs and respiratory ICWs of the AP-HP, Paris, France. The 18 participating centers (15 medical ICUs and affiliated ICWs, and three respiratory ICWs) are listed in Additional file 1: Table S1. All patients with severe Influenza virus infection consecutively admitted to the participating centers during the 2017-2018 influenza season were identified using the medical information system coding database (Programme de Médicalisation des Systèmes d'Information [PMSI]). The selection of adult stays (15 years and over) was performed on Diagnosis Related Groups in ICU/ICW, with the mention of "Influenza" in one of the coded diagnoses, using the International Classification of Diseases ICD-10 diagnosis codes (see Additional file 1). The diagnosis of Influenza virus infection was definite or probable, whether it was eventually microbiologically confirmed or not. Inter-institutional transfers within the AP-HP centers were considered in the patient care pathway by grouping together patient stays to obtain a database of unique patients.

Data recorded

For each selected case, baseline demographics and comorbidities, initial clinical presentation and vital signs, therapeutic management, ICU and hospital lengths of stay and vital status at discharge were extracted from the electronic health records (see Additional file 1).

Statistics

The characteristics of the population are described and compared according to their vital status at day 28 (D28). Qualitative variables are described by their frequencies and percentages of observed values, quantitative variables by their medians and interquartile ranges (IQR). Variables associated with 28-day mortality were identified using univariable Cox regression, with follow-up censored on D28. Hazard ratios (HR) are reported with their 95% confidence interval (CI). A sensitivity analysis was conducted to account for the multicenter design with the use of frailty models. Multivariable models were built to identify factors independently associated with 28-day mortality. Variables were included in the multivariable analysis from a practical perspective, when the information they provided was deemed clinically relevant and easily available on admission: two models were built from age > 65 years, comorbid conditions, abnormal core temperature (less than 35 °C or at least 40 °C), antiviral treatment on admission and a severity score [either the acute organ failure score (SOFA) or the CURB65 score]. Age > 65 years was removed from the model with the CURB65 score, as it already was a component of this score. When values were missing for components needed to compute the PSI or CURB65 scores [START_REF] Fine | A prediction rule to identify low-risk patients with community-acquired pneumonia[END_REF][START_REF] Lim | Defining community acquired pneumonia severity on presentation to hospital: an international derivation and validation study[END_REF], multiple imputation of these variables was used to compute these scores for all patients. Individual imputed patient's scores were averaged and rounded over 30 imputed datasets. Regression results relying on these scores were obtained by applying Rubin's rule on these datasets. All tests were two-tailed and p values < 0.05 were considered significant. Statistical analysis was conducted with R version 3.6.3 (R Core Team 2019; R foundation for statistical Computing, Vienna, Austria).

Ethical considerations

The EPIcuFLU_APHP research is a multicenter noninterventional data-based research using the care data collected during patients' stays at AP-HP. It was approved by the Scientific and Ethics Committee of the EDS AP-HP, which was authorized by the National Commission on Informatics and Liberty (CNIL) for such a noninterventional data-based research with no informed consent. There is no processing of indirectly identifiable data, or chaining with data from other sources, or longterm patient follow-up for this research. The access to data from different units and services of the AP-HP was the subject of the requesting investigator's fair information to data-producing professionals including the department heads or their representatives to ensure that they did not object to the use of patient's data they had taken care of.

Results

During the study period, 320 patients with probable or confirmed influenza infection were admitted to the ICU/ICW of the participating centers, a median of 3 [START_REF] Petrova | The evolution of seasonal influenza viruses[END_REF][START_REF] Nelson | The evolution of epidemic influenza[END_REF][START_REF] Harper | Seasonal influenza in adults and children-diagnosis, treatment, chemoprophylaxis, and institutional outbreak management: clinical practice guidelines of the Infectious Diseases Society of America[END_REF][START_REF] Delabre | Immunity against influenza A(H1N1) infections is determined by age at the time of initial strain circulation[END_REF][START_REF] Lemaitre | Mortality burden of the 2009 A/H1N1 influenza pandemic in France: comparison to seasonal influenza and the A/H3N2 pandemic[END_REF] days after symptoms onset, and 76 patients (24%) had already started taking oseltamivir before ICU/ICW admission (Table 1). The patients (188 men; 58.8%) were aged 63.2 [52. 3-73.4] years, had moderate overweight (body mass index > 30, n = 29; 9.1%) and often comorbid conditions, mainly congestive heart failure (n = 54; 16.9%), chronic renal disease (n = 53; 16.6%), and neoplastic disease (n = 42; 13.1%). Most patients had acute respiratory failure (n = 211; 65.9%) on ICU/ICW admission. The SAPSII score and SOFA score were 37 [28-55] and 5 [START_REF] Nelson | The evolution of epidemic influenza[END_REF][START_REF] Harper | Seasonal influenza in adults and children-diagnosis, treatment, chemoprophylaxis, and institutional outbreak management: clinical practice guidelines of the Infectious Diseases Society of America[END_REF][START_REF] Delabre | Immunity against influenza A(H1N1) infections is determined by age at the time of initial strain circulation[END_REF][START_REF] Lemaitre | Mortality burden of the 2009 A/H1N1 influenza pandemic in France: comparison to seasonal influenza and the A/H3N2 pandemic[END_REF][START_REF] Mercat | Severe H1N1 2009 influenza infection in adults: the French experience[END_REF][START_REF] Boddington | Developing a system to estimate the severity of influenza infection in England: findings from a hospital-based surveillance system between 2010/2011 and 2014/2015[END_REF][START_REF] Ercole | Modelling the impact of an influenza A/H1N1 pandemic on critical care demand from early pathogenicity data: the case for sentinel reporting[END_REF], respectively. About 13% and 3% of the data were missing for the variables needed to compute the PSI and CURB65 scores, respectively. Most patients were in the highest PSI risk classes (PSI IV-V: n = 262; 81.9%), and 146 patients (45.7%) had a CURB65 higher than 2 (Table 2). At least one factor targeted by the vaccination recommendations (Additional file 1: Table S2) was identified in 270 patients (84.4%), but only 52/245 patients (21.2%) in whom this information was reported had been vaccinated. Diagnostic PCR tests were mainly performed on nasopharyngeal swabs or nasopharyngeal aspirates (80%). Influenza infection was laboratory-confirmed in 196 patients (61.3%). The viruses A (H1N1) pdm09 and B/Yamagata were the main circulating viruses (Table 3, Fig. 1).

Altogether, vital support was required in up to 60% of critically ill patients during ICU/ICW stay, including mechanical ventilation (n = 182; 56.9%), vasopressors (n = 100; 31.2%), and renal replacement therapy (n = 39; 12.2%) (Table 4). Additional therapies were administered in the most severe patients, including oseltamivir (n = 278; 87%), steroids (n = 67; 20.9%), inhaled nitric oxide (NO, n = 17; 5.3%) and extracorporeal membrane oxygenation (ECMO, n = 6; 1.9%). Lengths of ICU and hospital stay were 8 [START_REF] Delabre | Immunity against influenza A(H1N1) infections is determined by age at the time of initial strain circulation[END_REF][START_REF] Lemaitre | Mortality burden of the 2009 A/H1N1 influenza pandemic in France: comparison to seasonal influenza and the A/H3N2 pandemic[END_REF][START_REF] Mercat | Severe H1N1 2009 influenza infection in adults: the French experience[END_REF][START_REF] Boddington | Developing a system to estimate the severity of influenza infection in England: findings from a hospital-based surveillance system between 2010/2011 and 2014/2015[END_REF][START_REF] Ercole | Modelling the impact of an influenza A/H1N1 pandemic on critical care demand from early pathogenicity data: the case for sentinel reporting[END_REF][START_REF]Equipes de surveillance de la grippe B-S Sibylle. Influenza activity in France, season 2017-2018[END_REF][START_REF] Fine | A prediction rule to identify low-risk patients with community-acquired pneumonia[END_REF][START_REF] Lim | Defining community acquired pneumonia severity on presentation to hospital: an international derivation and validation study[END_REF][START_REF] Gall | A new simplified acute physiology score based on a European-North American multicenter study[END_REF][START_REF] Vincent | The SOFA (Sepsis-related Organ Failure Assessment) score to describe organ dysfunction/failure. On behalf of the Working Group on sepsis-related problems of the European Society of Intensive Care Medicine[END_REF][START_REF] Pierre | Epidémiologie et de Santé Publique. Réseau Sentinelles[END_REF][START_REF] Adlhoch | The European Influenza Surveillance Network. Dominant influenza A(H3N2) and B/ Yamagata virus circulation in EU/EEA, 2016/17 and 2017/18 seasons, respectively[END_REF] days and 15 days, respectively.

At D28, 63 patients (19.7%) had died, 33 (10.3%) were still hospitalized in the ICU, while 224 (70.0%) had been discharged to conventional wards or long-term rehabilitation care units. The distribution of ICU admissions relative to the incidence of influenza-like illness in the Paris area is shown Fig. 1. Variables associated with 28-day-mortality are listed in Table 5, with hazard ratios for the risk of death 28 days after ICU/ICW admission. Two multivariable models were built, attempting at providing pragmatic and easily available information for clinical routine use. The first model identified four factors available during the first 24 h of ICU/ICW admission and independently associated with 28-day mortality: age > 65 years (HR: 1.79, 95% CI 1.02-3.16 p = 0.043), core temperature < 35 °C or ≥ 40 °C (HR: 3.06, 95% CI 1.73-5.42; p < 0.001), acute organ failure score (HR: 1.16, 95% CI 1.10-1.23 per 1-point increase; p < 0.001), and antiviral treatment on admission (HR: 0.45, 95% CI 0.24-0.85; p = 0.014). The scores dedicated to pneumonia (PSI and CURB65 scores) were not entered into that first model. The second model identified the CURB65 score > 2 (HR: 1.30, 95% CI 1.06-1.60; p = 0.014), core temperature < 35 °C or ≥ 40 °C (HR: 2.94, 95% CI 1.67-5.19; p < 0.001), and antiviral treatment on admission (HR: 0.34, 95% CI 0.19-0.61; p < 0.001) as being independently associated with 28-day mortality. Both these models were adjusted on the presence of comorbid conditions, which was not associated with mortality. The cumulative incidence of deaths over the 28-day followup is shown in Fig. 2, overall (A) and according to the CURB65 score (B), abnormal core temperature (C) and antiviral treatment on ICU admission (D). Frailty models accounting for the multicenter design yielded similar results (not reported).

Discussion

This non-interventional research was conducted during the flu season 2017-2018 in the adult ICUs/ICWs network of the AP-HP to assess the burden of the

Table 1 Characteristics of the population on ICU/ICW admission

Percentages were calculated after removing observations with missing data CGS Coma Glasgow scale, ICU intensive care unit, ICW intermediate care wards, MV mechanical ventilation, NSAIDs non-steroidal anti-inflammatory drugs a Obesity defined as Body Mass Index (kg/m 2 ) > 30 kg/m 2 b History of influenza vaccination was the variable with the largest amount of missing data; the reported percentages are 43/(257 -57) = 21.5%, 9/(63 -18) = 20.0% and 52/(320 -75) = 21.2% c Pregnancy, obesity (BMI > 30), 65-year old subjects and over, nursing home residency regardless of age, type 1 and 2 diabetes, chronic respiratory disease (chronic broncho-pulmonary diseases including asthma, broncho-pulmonary dysplasia and cystic fibrosis, chronic respiratory insufficiency), cardiac disease (congenital heart disease, heart failure, valvular disease, severe arrhythmia, coronary disease), neurological or muscle disease (stroke, severe forms of neurological and muscular disorders, para and tetraplegia with diaphragmatic involvement), renal disease (severe chronic renal insufficiency, nephrotic syndrome), immunosuppressive state (primary or acquired immune deficiency, except regular treatment with immunoglobulins, HIV infection and AIDS, solid transplantation), and others (hepatopathy, sickle cell disease, healthcare professionals) (see Additional file 1: Table S2) epidemic in critically ill patients, using the healthcare data warehouse of the Paris university hospitals network Early prognostic factors were mainly related to age, fever, the severity of the acute infection, as well as to the early administration of antiviral treatment. As expected, both the generic severity scores and scores dedicated to community-acquired pneumonia were strongly associated with 28-day mortality. The two final selected models highlight the importance of age and acute organ failure, two prognostic elements readily available at initial evaluation of patients. On the basis of the 2016-2017 flu season in France, we had estimated that approximately Table 5 Hazard ratios for the risk of death up to day 28 after admission to ICU/ICW * p < 0.05; **p < 0.001 a Pregnancy, obesity (body mass index > 30 kg/m 2 ), 65-year-old subjects and over, nursing home residency regardless of age, type 1 and 2 diabetes, chronic respiratory disease (chronic broncho-pulmonary diseases including asthma, broncho-pulmonary dysplasia and cystic fibrosis, chronic respiratory insufficiency), cardiac disease (congenital heart disease, heart failure, valvular disease, severe arrhythmia, coronary disease), neurological or muscle disease (stroke, severe forms of neurological and muscular disorders, para and tetraplegia with diaphragmatic involvement), renal disease (severe chronic renal insufficiency, nephrotic syndrome), immunosuppressive status (primary or acquired immune deficiency, except regular treatment with immunoglobulins, HIV infection and AIDS, solid transplantation), and others (liver disease, sickle cell disease, healthcare professionals). See Additional file 1: of cases, respectively. The initial severity was high, as demonstrated by the presence of an acute organ failure in two-thirds of the cases, involving mainly the lungs. The generic scores [START_REF] Gall | A new simplified acute physiology score based on a European-North American multicenter study[END_REF][START_REF] Vincent | The SOFA (Sepsis-related Organ Failure Assessment) score to describe organ dysfunction/failure. On behalf of the Working Group on sepsis-related problems of the European Society of Intensive Care Medicine[END_REF] as well as the scores dedicated to community-acquired pneumonia [START_REF] Fine | A prediction rule to identify low-risk patients with community-acquired pneumonia[END_REF][START_REF] Lim | Defining community acquired pneumonia severity on presentation to hospital: an international derivation and validation study[END_REF] were concordant with that severity, suggesting an overall probability of death ranging from 20 to 30%. The observed 28-day mortality was 19.7% (95% CI 15.5-24.5). These findings are in accordance with different sources of data collected in France for the 2017-2018 season [START_REF] Pierre | Epidémiologie et de Santé Publique. Réseau Sentinelles[END_REF]. These surveillance data showed that the epidemic started early, had a significant severity, and was exceptionally long in the context of insufficient vaccination coverage in France, and suboptimal vaccine efficacy. The atypical dynamics was related to the successive circulation of the A (H1N1) pdm09 and B/Yamagata viruses and may have contributed to those findings [START_REF] Pierre | Epidémiologie et de Santé Publique. Réseau Sentinelles[END_REF][START_REF] Adlhoch | The European Influenza Surveillance Network. Dominant influenza A(H3N2) and B/ Yamagata virus circulation in EU/EEA, 2016/17 and 2017/18 seasons, respectively[END_REF]. Antiviral therapy should be initiated as early as possible, the earlier initiation being more likely to provide benefit. However, randomized controlled trial data are not available to assess the impact of oseltamivir use among hospitalized patients with severe disease [START_REF] Hurt | Debate regarding oseltamivir use for seasonal and pandemic influenza[END_REF]. When replacing oseltamivir on ICU/ICW admission by oseltamivir started before ICU/ICW admission in the selected multivariable models, the hazard ratios for death within 28 days were still in favor of a lower mortality, yet the difference was not significant (HR: 0.65 [95% CI 0.34-1.24] and 0.62 [0.32-1.19] for the models with age > 65 years and with the CURB65 score, respectively). These results suggest that even though antiviral treatment is certainly beneficial for these patients when taken before ICU/ICW admission (which our study was not designed to evaluate), its benefit remains likely even when taken early on ICU/ICW admission. Thus, our findings confirm the high impact of the disease for the population at risk and strengthen the need for prevention, especially by promoting higher vaccine coverage among people at risk and compliance to control measures to limit the spread of the virus; our study also strongly suggest the potential benefit of early antiviral treatment for subjects at risk.

Limitations

The limitations of our clinical research project are related to the fact that it was one of the very first designed to use the Assistance Publique-Hôpitaux de Paris clinical data warehouse. Research on health data warehouses has become increasingly popular in the last few years and such databases may be considered as a major source of information for clinical research in a near future. This opportunity obviously relies on the possibility to effectively extract meaningful information from this massive amount of data that were collected without a predefined research question or even a research objective. This paradigm is far from the traditional clinical research settings that are mostly used to answer such questions, and there is still a lot to do to transform these databases into the powerful research tools they may become. In our study, some information were lacking, such as the diagnosis of early bacterial co-infection or the do-not-resuscitate order decisions that may have had an impact on the outcome of patients. The reproducibility and the validity of the results provided by the Assistance Publique-Hôpitaux de Paris clinical data warehouse are a major concern, and should be evaluated in different fields, for instance by comparing them with the results obtained from "goldstandard" dedicated studies relative to a specific research question. Since we started this study, a lot has been made to enhance the opportunity to use the Assistance Publique-Hôpitaux de Paris clinical data warehouse for research, and new challenges are constantly addressed so that this database may be part of the clinical researcher's toolbox within a few years. New healthcare data types are continuously integrated into the warehouse, which now contains data on ICD-10 diagnosis codes, medical procedures, laboratory results including microbiology, imaging, medication dispensing and medical documents regarding over 10 million patients in 39 hospitals.

To summarize, the main difficulties of this pilot work were linked to the fact that the healthcare data warehouse was not designed a priori for real-time monitoring. However, coding of diagnoses and medical procedures over time in critical care units and regular updates and extractions every 24 or 48 h are possible and could help achieve this goal. Data directly extracted from electronic medical records provide very useful and detailed clinical and care pathway information with more precision than those of the current regional surveillance network. These data should be available in a timely fashion with the aim of providing situational awareness regarding the most severe cases and thereby improving their detection. In the era of emerging infectious diseases and pandemics, the development of tools to monitor easily and in real time the progression of an epidemic, its severity and its impact on critical care organizations is a public healthcare imperative [START_REF] Richard | Interest of a simple on-line screening registry for measuring ICU burden related to an influenza pandemic[END_REF][START_REF] Anzic Influenza | Critical care services and the H1N1 (2009) influenza epidemic in Australia and New Zealand in 2010: the impact of the second winter epidemic[END_REF][START_REF] Gomersall | Expanding ICU facilities in an epidemic: recommendations based on experience from the SARS epidemic in Hong Kong and Singapore[END_REF][START_REF] Lum | Impact of pandemic (H1N1) 2009 influenza on critical care capacity in Victoria[END_REF]. Our pilot study is of particular relevance in the context of emerging respiratory viral pandemics, such as the COVID19 pandemic. It supports the usefulness of institutional tools to monitor the burden of the most severe cases in real time and to inform critical care services and health authorities, to adapt the healthcare system in a timely fashion by generating more ICU/ ICW capacities, sensitizing the emergency departments and finally contributing to the recommendations from health authorities. statistical analysis. MF and NL drafted the manuscript. All the authors critically revised the manuscript for important intellectual content. All the authors read and approved the final manuscript.

Fig. 1

 1 Fig. 1 ICU/ICW admissions throughout the 2017-2018 influenza season according to virus type. ICU/ICW weekly admissions in ICU/ICW, according to identified influenza virus types (bar plot, left axis) and incidence of influenza-like illness in Paris area (dashed line, right axis) [reference: https:// www. senti web. fr/]

Fig. 2

 2 Fig. 2 Cumulative incidence of deaths following ICU/ICW admission, overall (a), and according to the CURB65 score (b), core temperature (c) and antiviral treatment on admission (d).Panels illustrate the combination of the independant factors associated with 28-day mortality (second model), including the class of CURB65 (which may be either > 2 or ≤ 2; panel a), the core temperature on ICU admission (which may be either < 35 °C or ≥ 40 °C, or ≥ 35 °C and < 40 °C; panel c), and the antiviral treatment administered on ICU admission (panel d)

Table 2

 2 Severity of the disease during the first 24 hours of ICU/ICW admission

	Variable	All participants	Survivors	Non-survivors
		n = 320	n = 257	n = 63
	Acute organ failure
	Respiratory failure
	RR > 30/min or mechanical ventilation, n (%)	211 (65.9)	165 (64.0)	46 (73.8)
	PaO 2 /FiO 2 ratio (mmHg), median [IQR]	187 [90-336]	200 [98-350]	142 [73-260]
	ARDS, n (%)	120 (37.7)	82 (32.0)	38 (61.0)
	Shock
	SBP < 90, DBP < 60 or vasopressors, n (%)	116 (36.2)	81 (31.3)	35 (56.0)
	Renal failure
	Renal replacement therapy, n (%)	22 (6.9)	12 (4.7)	10 (15.9)
	Altered mental status (CGS < 14), n (%)	83 (28.2)	66 (27.8)	17 (29.8)
	Severity scores
	SAPS II (points), median [IQR]	37 [28-55]	34 [26-48]	59 [40-81]
	SOFA (points), median [IQR]	5 [2-8]	4 [2-7]	8 [4-14]
	Pneumonia Severity Index, median [IQR]	129 [100-160]	122 [96-153]	152 [131-186]
	Risk class, n (%)
	II	35 (7.9)	25 (9.7)	0 (0.0)
	III	33 (10.2)	31 (11.9)	2 (3.3)
	IV	105 (32.8)	92 (35.7)	13 (21.4)
	V	157 (49.1)	110 (42.8)	47 (75.2)
	CURB65, median [IQR]	2 [1-3]	2 [1-3]	3 [2-4]
	Risk class, n (%)
	0	24 (7.4)	21 (8.2)	3 (4.1)
	1	67 (21.0)	58 (22.6)	9 (14.5)
	2	83 (26.0)	72 (28.1)	11 (17.2)
	3	86 (26.7)	69 (26.8)	17 (26.3)
	4	41 (12.9)	26 (10.0)	15 (15.4)
	5	19 (6.1)	11 (4.3)	8 (8.4)

Components of the severity scores were imputed when missing: RR > 30/min or mechanical ventilation (n = 6), PaO 2 /FiO 2 ratio (n = 96), ARDS (n = 9), shock (n = 6), altered mental status (n = 26), SAPS II (n = 6). Patient counts for these scores were rounded over the 30 imputed datasets RR respiratory rate, PaO 2 partial pressure of oxygen, FiO 2 fraction of inspired oxygen, ARDS acute respiratory distress syndrome, SBP/DBP systolic/diastolic blood pressure, CGS Glasgow coma scale, SAPS simplified acute physiology score, SOFA sequential organ failure assessment

Table 3

 3 Microbiological diagnosis

	Variable	All participants	Survivors	Non-survivors
		n = 320	n = 257	n = 63
	Virus sampling
	Nasopharyngeal swab	249 (77.8)	203 (79.0)	46 (73.0)
	Nasopharyngeal aspirate	10 (3.1)	7 (2.7)	3 (4.8)
	Sputum	14 (4.4)	10 (3.9)	4 (6.3)
	Bronchial aspirate	31 (9.7)	23 (8.9)	8 (12.7)
	Influenza virus type and subtype identified in patients with at least one positive test	n = 196	n = 155	n = 41
	Virus A a	78 (39.8)	59 (38.1)	19 (46.3)
	A/H3N2 a	10 (5.1)	6 (3.9)	4 (9.8)
	A/H1N1pdm a	68 (34.7)	53 (32.2)	15 (36.6)
	Virus B a	120 (61.2)	96 (61.9)	24 (58.6)

a Positive tests are reported as counts and percentages in patients with at least one positive test (EDS AP-HP). Our findings highlight the high burden and severity of influenza on critical care services, involving patients with low vaccination coverage, requiring life support in 60% of cases, and having prolonged length of stay and a high (20%) mortality rate.

Table 4

 4 Treatments administered in ICU and outcomes NO nitric oxide, ECMO extracorporeal membrane oxygenation a Invasive mechanical ventilation with endotracheal intubation, non-invasive positive pressure ventilation and high flow nasal oxygen therapies were included in the category mechanical ventilation b Steroids included hydrocortisone and methylprednisolone c Oseltamivir was either continued (n = 76) or introduced within the first 24 h of ICU/ICW referral (n = 202)

	Variable	All participants	Survivors	Non-survivors
		n = 320	n = 257	n = 63
	Organ support and additional therapies, n (%)
	Mechanical ventilation a	182 (56.9)	136 (52.9)	46 (73.0)
	Invasive mechanical ventilation	116 (36.2)	77 (30.0)	39 (61.9)
	Neuromuscular blocking agents	17 (5.3)	13 (5.1)	4 (6.3)
	Prone positioning	35 (10.9)	23 (8.9)	12 (19.0)
	NO	17 (5.3)	13 (5.1)	4 (6.3)
	ECMO	6 (1.9)	3 (1.2)	3 (4.8)
	Vasopressors	100 (31.2)	62 (24.1)	38 (60.3)
	Renal replacement therapy	39 (12.2)	22 (8.6)	17 (27)
	Steroids b	67 (20.9)	52 (20.2)	15 (23.8)
	Oseltamivir c	278 (86.9)	232 (90.3)	46 (73.0)
	Length of stay (days), median [IQR]
	ICU	8 [4-15]	8 [4-15]	7 [2-15.5]
	Hospital	15 [8-28]	16 [9-31]	10 [3.5-18.5]

Table S2

 S2

	Variable	Univariable analysis		Multivariable analysis
		HR (95% CI)	p	First model	Second model
				HR (95% CI)	HR (95% CI)
	Age (per 10 years)	1.42 (1.18-1.71)	0.0002	-	-
	Age > 65 years	1.70 (1.03-2.81)	0.039	1.79 (1.02-3.16)*	-
	Sex (male)	1.20 (0.72-2.00)	0.49	-	-
	Obesity	1.11 (0.48-2.58)	0.81	-	-
	Comorbid conditions	1.33 (0.81-2.18)	0.26	1.16 (0.69-1.93)	1.17 (0.70-1.97)
	Neoplastic disease	1.50 (0.80-2.82)	0.20	-	-
	Congestive heart failure	1.00 (0.52-1.92)	0.99	-	-
	Cerebrovascular disease	1.78 (0.72-4.45)	0.21	-	-
	Renal disease	1.36 (0.75-2.46)	0.31	-	-
	Liver disease	6.64 (2.05-21.54)	0.0016	-	-
	Status of influenza vaccination	1.03 (0.48-2.21)	0.94	-	-
	At least one factor targeted by the vaccination a	1.87 (0.75-4.67)	0.18	-	-
	Time to referral from symptoms onset (per day)	1.00 (0.97-1.03)	0.99	-	-
	Treatments before admission	0.89 (0.52-1.52)	0.67
	Antimicrobial drugs	0.81 (0.49-1.34)	0.41	-	-
	Antibiotics	0.85 (0.51-1.39)	0.51	-	-
	Antiviral therapy (oseltamivir)	0.70 (0.37-1.31)	0.26	-	-
	Anti-inflammatory drugs	0.89 (0.44-1.8)	0.75	-	-
	Steroids	0.74 (0.32-1.7)	0.47	-	-
	NSAIDs	1.62 (0.59-4.41)	0.35	-	-
	Physical exam during the first 24 h of ICU/ICW admission
	Pulse ≥ 125 per minute Temperature < 35 °C or ≥ 40 °C Respiratory rate ≥ 30 per minute or MV Systolic blood pressure < 90 mm Hg or diastolic blood pres-	1.21 (0.71-2.06) 3.37 (1.96-5.81) 1.49 (0.85-2.62) 2.4 (1.44-3.99)	0.48 0.00002 0.17 0.0007	-3.06 (1.73-5.42)** --	-2.94 (1.67-5.19)** --
	sure < 60 mm Hg or vasopressors
	Renal replacement therapy	2.6 (1.32-5.11)	0.0056	-	-
	Altered mental status (CGS < 14)	1.05 (0.6-1.84)	0.87	-	-
	Initial severity (scores)
	SAPS2 (points) 300 patients would be hospitalized in the participating 1.04 (1.03-1.05)	< 10e-5	-	-
	SOFA (points) adult ICUs over a similar 2017-2018 season, i.e., with a 1.15 (1.09-1.21)	< 10e-5	1.16 (1.10-1.23)**	-
	Pneumonia Severity Index, Class IV-V "severe" profile, and that the expected number of deaths 7.30 (1.04-51.31)	0.046	-	-
	CURB65 > 2 at 28 days would range between 43 (moderate epidemic) 2.20 (1.30-3.75)	0.0035	-	1.30 (1.06-1.60)*
	Antiviral therapy and 86 (severe epidemic). Altogether, we recorded 320
	Oseltamivir (ongoing or introduced after ICU admission) ICU stays involving 65-year or older and fragile patients 0.32 (0.18-0.56)	0.00017	0.45 (0.24-0.85)*	0.34 (0.19-0.61)**
	with comorbid conditions, and 73 (22.8%) in-hospital
	deaths. At least one factor targeted by the vaccination
	recommendations was identified in most of them, but
	only 20% had received vaccination, reinforcing the fact
	that efforts to foster preventive strategies including vac-
	cination are needed. Most patients were directly referred
	to the ICUs from out-of-hospital emergency services or
	EDs, about 3 days after symptoms onset. The patients
	were receiving antimicrobial therapies (antibiotics and/
	or oseltamivir) and/or anti-inflammatory drugs (steroids
	and/or NSAIDs) before ICU referral in 61.3% and 16.6%

Acknowledgements

The authors thank EPIcuFLU_APHP group.

Funding

The authors have no funding for the research to declare.

Availability of data and materials

The dataset used and analyzed during the current study are available from the AP-HP Clinical Data Warehouse on reasonable request.

Abbreviations

Supplementary Information

The online version contains supplementary material available at https:// doi. org/ 10. 1186/ s13613-021-00884-8.

Additional file 1.

Selection of the population and Data recorded. Table S1. Participating centres. Table S2. Factor(s) targeted by the vaccination in the studied population.

Authors' contributions

NL and FC had full access to all of the data in the study and take responsibility for the integrity of the data and the accuracy of the data analysis. MF, NL and FC contributed to study concept and design. All the authors helped in acquisition of data. MF, NL and FC analyzed and interpreted the data. NL made the

Declarations Ethics approval and consent to participate

The EPIcuFLU_APHP was approved by the Scientific and Ethics Committee of the EDS AP-HP, which was authorized by the National Commission on Informatics and Liberty (CNIL) for such a non-interventional data-based research with no informed consent.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests to declare.

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.