Gaëtan Heidsieck

Daniel De Oliveira

Esther Pacitti

Christophe Pradal

François Tardieu

Patrick Valduriez

Distributed Caching of Scientific Workflows in Multisite Cloud

Keywords: Multisite cloud, Distributed Caching, Scientific Workflow, Workflow System, Workflow Scheduling

Many scientific experiments are performed using scientific workflows, which are becoming more and more data-intensive. We consider the efficient execution of such workflows in the cloud, leveraging the heterogeneous resources available at multiple cloud sites (geo-distributed data centers). Since it is common for workflow users to reuse code or data from other workflows, a promising approach for efficient workflow execution is to cache intermediate data in order to avoid re-executing entire workflows. In this paper, we propose a solution for distributed caching of scientific workflows in a multisite cloud. We implemented our solution in the OpenAlea workflow system, together with cache-aware distributed scheduling algorithms. Our experimental evaluation on a three-site cloud with a data-intensive application in plant phenotyping shows that our solution can yield major performance gains, reducing total time up to 42% with 60% of same input data for each new execution.

Introduction

In many scientific domains, e.g., bio-science [START_REF] Kelling | Data-intensive science: a new paradigm for biodiversity studies[END_REF], complex numerical experiments typically require many processing or analysis steps over huge datasets. They can be represented as scientific workflows, or workflows, for short, which facilitate the modeling, management and execution of computational activities linked by data dependencies. As the size of the data processed and the complexity of the computation keep increasing, these workflows become data-intensive [START_REF] Kelling | Data-intensive science: a new paradigm for biodiversity studies[END_REF], thus requiring high-performance computing resources.

The cloud is a convenient infrastructure for handling workflows, as it allows leasing resources at a very large scale and relatively low cost. In this paper, we consider the execution of a large workflow in a multisite cloud, i.e., a cloud with geo-distributed cloud data centers (sites). Note that a multisite option is now well supported by all popular public clouds, e.g., Microsoft Azure, Amazon EC2, and Google Cloud, which provide the capability of using multiple sites with a single cloud account, thus avoiding the burden of using multiple accounts. The main reasons for using multiple cloud sites for data-intensive workflows is that they often exceed the capabilities of a single site, either because the site imposes usage limits for fairness and security, or simply because the datasets are too large. In scientific applications, there can be much heterogeneity in the storage and computing capabilities of the different sites, e.g., on premise servers, HPC platforms from research organizations or federated cloud sites at the national level [START_REF] Crago | Heterogeneous cloud computing[END_REF]. As an example in plant phenotyping, greenhouse platforms generate terabytes of raw data from plants, which are typically stored at data centers geographically close to the greenhouse to minimize data transfers. However, the computation power of those data centers may be limited and fail to scale when the analyses become more complex, such as in plant modeling or 3D reconstruction. Other computation sites are then required.

Most Scientific Workflow Management Systems (workflow systems) can execute workflows in the cloud [START_REF] De Oliveira | Towards a taxonomy for cloud computing from an e-science perspective[END_REF]. Some examples are Swift/T, Pegasus, SciCumulus, Kepler and OpenAlea [START_REF] Liu | A survey of data-intensive scientific workflow management[END_REF]. Our work is based on OpenAlea [START_REF] Pradal | Openalea: scientific workflows combining data analysis and simulation[END_REF], which is widely used in plant science for simulation and analysis. Most existing systems use naive or user-based approaches to distribute the tasks across sites. The problem of scheduling a workflow execution over a multisite cloud has started to be addressed in [START_REF] Maheshwari | Improving multisite workflow performance using model-based scheduling[END_REF], using performance models to predict the execution time on different resources. In [START_REF] Liu | Multi-objective scheduling of scientific workflows in multisite clouds[END_REF], we proposed a solution based on multi-objective scheduling and a single site virtual machine provisioning approach, assuming homogeneous sites, as in public cloud.

Since it is common for workflow users to reuse code or data from other workflows [START_REF] Garijo | Common motifs in scientific workflows: An empirical analysis[END_REF], a promising approach for efficient workflow execution is to cache intermediate data in order to avoid re-executing entire workflows. Furthermore, a user may need to re-execute a workflow many times with different sets of parameters and input data depending on the previous results generated. Fragments of the workflow, i.e. a subset of the workflow activities and dependencies, can often be reused. Another important benefit of caching intermediate data is to make it easy to share with other research teams, thus fostering new analyses at low cost.

Caching has been supported by some workflow systems, e.g., Kepler, Vis-Trails and OpenAlea. Kepler [START_REF] Altintas | Provenance collection support in the kepler scientific workflow system[END_REF] provides a persistent cache on the cloud, but at a single site, and does not support multisite. VisTrails [START_REF] Callahan | Vistrails: visualization meets data management[END_REF] provides a persistent cache, but only for local execution on a personal desktop. In [START_REF] Heidsieck | Adaptive caching for data-intensive scientific workflows in the cloud[END_REF], we proposed an adaptive caching method for OpenAlea that automatically determines the most suited intermediate data to cache, taking into account workflow fragments, but only in the case of a single cloud site. Another interesting single site method, also exploiting workflow fragments, is to compute the ratio between re-computation cost and storage cost to determine what intermediate data should be stored [START_REF] Yuan | A highly practical approach toward achieving minimum data sets storage cost in the cloud[END_REF]. All these methods are single site (centralized). The only distributed caching method for workflow execution in a multisite cloud we are aware of is restricted to hot metadata (frequently accessed metadata) [START_REF] Liu | Efficient scheduling of scientific workflows using hot metadata in a multisite cloud[END_REF], ignoring intermediate data.

Caching data in a multisite cloud with heterogeneous sites is much more complex. In addition to the trade-off between re-computation and storage cost at single sites, there is the problem of site selection for placing cached data. The problem is more difficult than data allocation in distributed databases [START_REF] Özsu | Principles of Distributed Database Systems[END_REF], which deals only with well-defined base data, not intermediate data produced by tasks. Furthermore, the scheduling of workflow executions must be cacheaware, i.e., exploit the knowledge of cached data to decide between reusing and transferring cached data versus re-executing the workflow fragments.

In this paper, we propose a distributed solution for caching of scientific workflows in a multisite cloud. Based on a distributed and parallel architecture [START_REF] Özsu | Principles of Distributed Database Systems[END_REF], composed of heterogeneous sites (including on premise servers and shared-nothing clusters), we propose algorithms for adaptive caching, cache site selection and dynamic workflow scheduling. We implemented our caching solution in OpenAlea, together with a multisite scheduling algorithm. Based on a real data-intensive application in plant phenotyping, we provide an extensive experimental evaluation using a cloud with three heterogeneous sites.

This paper is organized as follows. Section 2 presents our real use case in plant phenotyping. Section 3 introduces our workflow system architecture in multisite cloud. Section 4 describes our cache management solution. Section 5 gives our experimental evaluation. Finally, Section 6 concludes.

Use Case in Plant Phenotyping

In this section, we introduce a real use case in plant phenotyping that will serve as motivation for the work and basis for the experimental evaluation. In the last decade, high-throughput phenotyping platforms have emerged to allow for the acquisition of quantitative data on thousands of plants in well-controlled environmental conditions. For instance, the seven facilities of the French Phenome project 5 produce each year 200 Terabytes of data, which are various (images, environmental conditions and sensor outputs), multiscale and originate from different sites. Analyzing such massive datasets is an open, yet important, problem for biologists [START_REF] Tardieu | Plant phenomics, from sensors to knowledge[END_REF].

The Phenomenal workflow ( [START_REF] Artzet | Openalea.phenomenal: A workflow for plant phenotyping[END_REF]) has been developed in OpenAlea to analyze and reconstruct the geometry and topology of thousands of plants through time in various conditions. It is composed of nine fragments such as image binarization, 3D volume reconstruction, organ segmentation or intercepted light simulation. Different users can conduct different biological analyses by reusing some workflow fragments on the same dataset to test different hypotheses [START_REF] Heidsieck | Adaptive caching for data-intensive scientific workflows in the cloud[END_REF]. To save both time and resources, they want to reuse the intermediate results that have already been computed rather than recompute them from scratch.

The raw data comes from the Phenoarch platform, which has a capacity of 1,680 plants within a controlled environment (e.g., temperature, humidity, irrigation) and automatic imaging through time. The total size of the raw image dataset for one experiment is 11 Terabytes. To limit data movement, the raw data is stored at a server near to the experimental platform, with both data storage and computing resources. However, these computing resources are not enough to process a full experiment in a relatively short time. Thus, scientists who need to do a full experiment will execute the Phenomenal workflow at a more powerful site by transferring the raw data for each new analysis.

In this Phenomenal use case, the cloud is composed of heterogeneous sites, with both on premise servers close to the experimental platform and other more powerful cloud sites. The on premise server has high storage capacity and hosts the raw data. Other sites are used to computational intensive executions, with high-performance computing resources. On premise servers are used locally to execute some Phenomenal fragments that do not require powerful resources. In this case, one has to choose between transferring the raw data or some intermediate data to a powerful site or re-executing some fragments locally before transferring intermediate data. The trade-off between data re-computation and data transfer is complex in a multisite cloud with much heterogeneity. In particular, one needs to pay attention to cached data placement, so as to avoid bottlenecks on the most used intermediate data.

Multisite Cloud Workflow System Architecture

In this section, we present our workflow system architecture that integrates caching and reuse of intermediate data in a multisite cloud. We motivate our design decisions and describe our architecture in terms of nodes and components (see Figure 1), which are involved in the processing of workflows.

Our architecture capitalizes on the latest advances in distributed and parallel data management to offer performance and scalability [START_REF] Özsu | Principles of Distributed Database Systems[END_REF]. We consider a distributed cloud architecture with on premise servers, where raw data is produced, e.g., by a phenotyping experimental platform in our use case, and remote sites, where the workflow is executed. The remote sites (data centers) are sharednothing clusters, i.e., clusters of server machines, each with processor, memory and disk. We adopt shared-nothing as it is the most scalable and cost-effective architecture for big data analysis.

In the cloud, metadata management has a critical impact on the efficiency of workflow scheduling as it provides a global view of data location, e.g., at which nodes some raw data is stored, and enables task tracking during execution [START_REF] Liu | Efficient scheduling of scientific workflows using hot metadata in a multisite cloud[END_REF]. We organize the metadata in three repositories: catalog, provenance database and cache index. The catalog contains all information about users (access rights, etc.), raw data location and workflows (code libraries, application code). The provenance database captures all information about workflow execution. The cache index contains information about tasks and cache data produced, as well as the location of files that store the cache data. Thus, the cache index itself is small (only file references) and the cached data can be managed using the underlying file system. A good solution for implementing these metadata repositories is a key-value store, such as Cassandra6 , which provides efficient key-based access, scalability and fault-tolerance through replication in a shared-nothing cluster.

The raw data (files) are initially produced and stored at some cloud sites, e.g., in our use case, at the phenotyping platform. During workflow execution, the intermediate data is generated and consumed at one site's node in memory. It gets written to disk when it must be transferred to another node (potentially at the same site), or when explicitly added to the cache. The cached data (files) can later be replicated at other sites to minimize data transfers.

We extend the workflow system architecture proposed in [START_REF] Liu | A survey of data-intensive scientific workflow management[END_REF] for single site. It is composed of six modules: workflow manager, global scheduler, local scheduler, task manager, data manager and metadata manager, to support both execution and intermediate data caching in a multisite cloud. The workflow manager provides a user interface for workflow definition and processing. Before workflow execution, the user selects a number of virtual machines (VMs), given a set of possible instance formats, i.e., the technical characteristics of the VMs, deployed on each site's nodes. When a workflow execution is started, the workflow manager simplifies the workflow by removing some workflow fragments and partitions depending on the raw input data and the cached data (see Section 4). The global scheduler uses the metadata (catalog, provenance database, and cache index) to schedule the workflow fragments of the simplified workflow. The VMs on each site are then initialized, i.e., the programs required for the execution of the tasks are installed and all parameters are configured. The local scheduler schedules the workflow fragments received on its VMs. The data manager module handles data transfers between sites during execution (for both newly generated intermediate data and cached data) and manages cache storage and replication. At a single site, data storage is distributed between nodes. Finally, the task manager (on each VM) manages the execution of fragments on the VMs at each site. It exploits the provenance metadata to decide whether or not the task's output data should be placed in the cache, based on the cache provisioning algorithm described in Section 4. Local scheduling and execution can be performed as in [START_REF] Heidsieck | Adaptive caching for data-intensive scientific workflows in the cloud[END_REF].

Figure 1 shows how these components are involved in workflow processing, using the traditional master-worker model. In this architecture, we consider two types of cloud sites, i.e., coordinator and participant. The relationship between the site is also based on the master-worker model, the coordinator site, managing the participant sites. The workflow manager and the global scheduler modules are implemented on the coordinator site. The remaining modules are implemented on all sites.

At each site, there are three kinds of nodes: master, compute and data nodes, which are mapped to cluster nodes at configuration time, e.g. using a cluster manager like Yarn ( http://hadoop.apache.org). There is one active master node per site. There is also a standby node to deal with master node failure. The master nodes are the only ones to communicate across sites. The local scheduler and metadata management modules are implemented on the master node, which manages communication, metadata and scheduling. The master nodes are responsible for transferring data between sites during execution.

Multisite Cache-aware Workflow Execution

In this section, we present in more details how the global scheduler performs multisite cache-aware workflow execution. In particular, the global scheduler must decide which data to cache (cache data selection) and where (cache site selection), and where to execute workflow fragments (execution site selection). Since these decisions are not independent, we propose a cost function to make a global decision, based on the cost components for individual decisions. We start by giving an overview of distributed workflow execution. Then, we present the methods and cost functions for cache data selection, cache site selection and execution site selection. Finally, we introduce our cost function for the global decision.

Distributed Workflow Execution Overview

We consider a multisite cloud with a set of sites S={s 1 , ..., s n }. A workflow W (A, D) is a a directed acyclic graph (DAG) of computational activities A and their data dependencies D. A task t is the instantiation of an activity during execution with specific associated input data. A fragment f of an instantiated workflow is a subset of tasks and their dependencies.

The execution of a workflow W (A, D) in S starts at a coordinator site s c and proceeds in three main steps:

1. The global scheduler at s c simplifies and partitions the workflow into fragments. Simplification uses metadata to decide whether a task can be replaced by corresponding cached data references. Partitioning uses the dependencies in D to produce fragments.

2. For each fragment, the global scheduler at s c computes a cost function to make a global decision on which data to cache where, and on which site to execute. Then, it triggers fragment execution and cache placement at the selected sites. 3. At each selected site, the local scheduler performs the execution of its received fragments using its task manager (to execute tasks) and data manager (to transfer the required input data). It also applies the decision of the global scheduler on storing new intermediate data into the cache.

We introduce basic cost functions to reflect data transfer and distributed execution. The time to transfer some data d from site s i to site s j , noted T tr (d, s i , s j ), is defined by

T tr (d, s i , s j ) = Size(d) T rRate(s i , s j ) (1)
where T rRate(s i , s j ) is the transfer rate between s i and s j .

The time to transfer input and cached data, In(f ) and Cached(f ) respectively, to execute a fragment f at site s i is T input (f, s i ):

T input (f, s i ) = S sj (T tr (In(f ), s j , s i ) + T tr (Cached(f ), s j , s i )) (2)
The time to compute a fragment f at site s, noted T compute (f, s), can be estimated using Amdahl's law [START_REF] Zhang | Scheduling of scientific workflow in non-dedicated heterogeneous multicluster platform[END_REF]:

T compute (f, s) = ( α n + (1 -α)) * W (f ) CP U perf (s) (3)
where W (f ) is the workload for the execution of f , CP U perf (s) is the average computing performance of the CPUs at site s and n is the number of CPUs at site s. We suppose that the local scheduler may parallelize task executions. Therefore, α represents the percentage of the workload that can be executed in parallel.

The expected waiting time to be able to execute a fragment at site s is noted T wait (s), which is the minimum expected time for s to finish executing the fragments in its queue.

The time to transfer the intermediate data generated by fragment f at site s i to site s j , noted T write (Output(f ), s i , s j ), is defined by:

T write (Output(f ), s i , s j ) = T tr (Output(f ), s i , s j ) (4)
where Output(f ) is the data generated by the execution of f .

Cache Data Selection

To determine what new intermediate data to cache, we consider two different methods: greedy and adaptive. Greedy data selection simply adds all new data to the cache. Adaptive data selection extends our method proposed in [START_REF] Heidsieck | Adaptive caching for data-intensive scientific workflows in the cloud[END_REF] to multisite cloud. It achieves a good trade-off between the cost saved by reusing cached data and the cost incurred to feed the cache.

To determine if it is worth adding some intermediate data Output(f ) at site s j , we consider the trade-off between the cost of adding this data to the cache and the potential benefit if this data was reused. The cost of adding the data to site s j is the time to transfer the data from the site where it was generated. The potential benefit is the time saved from loading the data from s j to the site of computation instead of re-executing the fragment. We model this trade-off with the ratio between the cost and benefit of the cache, noted p(f, s i , s j ), which can be computed from equations 2, 3 and 4,

p(f, s i , s j ) = T write (Output(f ), s i , s j ) T input (f, s i ) + T compute (f, s i ) -T tr (Output(f ), s j , s i ) (5)
In the case of multiple users, the probability that Output(f ) will be reused or the number of times fragment f will be re-executed is not known when the workflow is executed. Thus, we introduce a threshold T hreshold (computed by the user) as the limit value to decide whether a fragment output will be added to the cache. The decision on whether Output(f ) generated at site s i is stored at site s j can be expressed by

i,j = 1, if p(f, s i , s j ) < T hreshold. 0, otherwise. (6)

Cache Site Selection

Cache site selection must take into account the data transfer cost and the heterogeneity of computing and storage resources. We propose two methods to balance either storage load (bStorage) or computation load (bCompute) between sites. The bStorage method allows preventing bottlenecks when loading cached data.

To assess this method at any site s, we use a load indicator, noted L bStorage (s), which represents the relative storage load as the ratio between the storage used for the cached data (Storage used (s)) and the total storage (Storage total (s)).

L bStorage (s) = Storage used (s) Storage total (s) (7)
The bCompute method balances the cached data between the most powerful sites, i.e., with more CPUs, to prevent computing bottlenecks during execution. Using the knowledge on the sites' computing resources and usage, we use a load indicator for each site s, noted L bCompute (s), based on CPUs idleness (CP U idle (s)) versus total CPU capacity (CP U total (s)).

L bCompute (s) = 1 -CP U idle (s) CP U total (s) (8)
The load of a site s, depending on the method used, is represented by L(s), ranging between 0 (empty load) and 1 (full). Given a fragment f executed at site s i , and a set of sites s j with enough storage for Output(f ), the best site s * to add Output(f ) to its cache can be obtained using Equation 1(to include transfer time) and Equation 6(to consider multiple users),

s * (f ) si = argmax sj ( i,j * (1 -L(s j )) T write (Output(f ), s i , s j ) ) (9)

Execution Site Selection

To select an execution site s for a fragment f , we need to estimate the execution time for f as well as the time to feed the cache with the result of f . The execution time f at site s (T execute (f, s)) is the sum of the time to transfer input and cached data to s, the time to get computing resources and the time to compute the fragment. It is obtained using Equations 2 and 3.

T execute (f, s) = T input (f, s) + T compute (f, s) + T wait (s) ( 10
)
Given a fragment f executed at site s i and its intermediate data Output(f ), the time to write Output(f ) to the cache (T f eed cache (f, s i , s j )) can be defined as:

T f eed cache (f, s i , s j , i,j ) = i,j * T write (Output(f ), s i , s j ) ( 11
)
where s j is given by Equation 9.

Global Decision

At Step 2 of workflow execution, for each fragment f , the global scheduler must decide on the best combination of individual decisions regarding cache data, cache site, and execution site. These individual decisions depend on each other. The decision on cache data depends on the site where the data is generated and the site where it will be stored. The decision on cache site depends on the site where the data is generated and the decision of whether or not the data will be cached. Finally, the decision on execution site depends on what data will be added to the cache and at which site. Using Equations 10 and 11, we can estimate the total time (T total ) for executing a fragment f at site s i and adding its intermediate data to the cache at another site s j :

T total (f, s i , s j , i,j ) = T execute (f, s i ) + T f eed cache (f, s i , s j , i,j ) (12)
Then, the global decision for cache data ( (f )), cache site (s * cache ) and execution site (s * exec ) is based on minimizing the following equation for the n 2 pairs of sites s i and s j

(s * exec , s * cache , (f )) = argmin si,sj (T total (f, s i , s j , i,j )) ( 13
)
This decision is done by the coordinator site at before each fragment execution and only takes into account the cloud site's status at that time. Note that s * exec , s * cache can be the coordinator site and can be the same site.

Experimental Evaluation

In this section, we first present our experimental setup, which features a heterogeneous multisite cloud with multiple users who re-execute part of the workflow. Then, we compare the performance of our multisite cache scheduling method against two baseline methods. We end the section with concluding remarks.

Experimental Setup

Our experimental setup includes a multisite cloud, with three sites in France, a workflow implementation and an experimental dataset. Site 1 in Montpellier is a server close to the Phenoarch phenotyping platform. It has the smallest number of CPUs and largest amount of storage among the sites. The raw data is stored at this site. Site 2 is the coordinator site, located in Lille. Site 3, located in Lyon, has the largest number of CPUs and the smallest amount of storage.

To model site heterogeneity in terms of storage and CPU resources, we use heterogeneity factor H in three configurations: H = 0, H = 0.3 and H = 0.7. For the three sites altogether, the total number of CPUs is 96 and the total storage on disk for intermediate data is 180 GB (The raw data is stored on an additional node at Site 1). On each site, several nodes are instantiated for the executions, they have a determined number of CPUs from 1, 2, 4, 8 or 16 CPUs. The available disk size for each node is limited by implementation. With H = 0 (homogeneous configuration), each site has 32 CPUs (two 16 CPUs nodes) and 60 GB (30 GB each). With H = 0.3, we have 22 CPUs and 83 GB for Site 1, 30 CPUs and 57 GB for Site 2 and 44 CPUs and 40 GB for Site 3. With H = 0.7 (most heterogeneous configuration), we have 6 CPUs and 135 GB for Site 1, 23 CPUs and 35 GB for Site 2 and 67 CPUs and 10 GB for Site 3.

The input dataset for the Phenomenal workflow is produced by the Phenoarch platform (see Section 2). Each execution of the workflow is performed on a subset of the input dataset, i.e. 200 GB of raw data, which represents the execution of 15,000 tasks. For each user, 60% of the raw data is reused from previous executions. Thus each execution requires only 40% of new raw data. For the first execution, no data is available in the cache.

We implemented our cache-aware scheduling method, which we call cacheA, in OpenAlea and deployed it at each site using the Conda multi-OS package manager. The metadata distributed database is implemented using Cassandra. Communication between the sites is done using the protocol library ZeroMQ. Data transfer between sites is done through SSH. We have also implemented two baseline methods, Sgreedy and Agreedy, based on the SiteGreedy and Act-Greedy methods described in [START_REF] Liu | Multi-objective scheduling of scientific workflows in multisite clouds[END_REF], respectively. The Sgreedy method extends SiteGreedy, which schedules each workflow fragment at a site that is available for execution, with our cache data/site selection methods. Similarly, the Agreedy method extends ActGreedy, which schedules each workflow fragment at a site that minimizes a cost function based on execution time and input data transfer time, with our cache data/site selection methods. These baseline methods perform execution site selection followed by cache data/site selection while CacheA makes a global decision.

Experiments

We compare CacheA with the two baseline methods in terms of execution time and amount of data transferred. We define total time as execution time plus transfer time. In experiment 1, we consider a workflow execution with caching or without. In Experiment 2, we consider multiple users who execute the same workflow on similar input data, where 60% of the data is the same. In Experiment 3, we consider different heterogeneous configurations for one workflow execution.

Experiment 1: with caching. In this basic experiment, we compare two workflow executions: with caching, using CacheA and bStorage; and without caching, using ActGreedy. We consider one re-execution of the workflow on different input datasets, from 0% to 60% of same reused data.

CacheA outperforms ActGreedy from 20% of reused data. Below 20%, the overhead of caching outweighs its benefit. For instance, with no reuse (0%), the total time with CacheA is 16% higher than with ActGreedy. But with 30%, it is 11% lower, and with 60%, it is 42% lower.

Experiment 2: multiple users. Figure 2 shows the total time of the workflow for the three scheduling methods, four users, H = 0.7 and our two cache site selection methods: (a) bStorage, and (b) bCompute.

Let us first analyze the results in Figure 2.a (bStorage method). For the first user execution, CacheA outperforms Sgreedy in terms of execution time by 8% and in terms of data and intermediate data transfer times by 51% and 63%, respectively. The reason Sgreedy is slower is that it schedules some computeintensive fragments at Site 1, which has the lowest computing resources. Furthermore, it does not consider data placement and transfer time when scheduling fragments.

Again for the first user execution, CacheA outperforms Agreedy in terms of total time by 24%, when considering data transfer time to the cache. However, CacheA execution time is a bit slower (by 9%). The reason that Agreedy is slower in terms of total time is that it does not take into account the placement of the cached data, which leads to larger amounts (by 67%) of cache data to transfer. For other users' executions (when cached data exists), CacheA outperforms Sgreedy in terms of execution time by 29%, and for the fourth user execution, by 20%. This is because CacheA better selects the cache site in order to reduce the execution time of the future re-executions. In addition, CacheA balances the cached data and computations. It outperforms Sgreedy and Agreedy in terms of intermediate data transfer times (by 59% and 15%, respectively), and cache data transfer times (by 82% and 74%, respectively). Overall, CacheA outperforms Sgreedy and Agreedy in terms of total times by 61% and 43%, respectively. The workflow fragments are not necessarily scheduled to the site with shortest execution time, but to the site that minimizes overall total time. Considering the multiuser perspective, CacheA outperforms baseline methods, reducing the total time for each new user (up to 6% faster for the fourth user compared to the second).

Let us now consider Figure 2.b (bCompute method). For the first user execution, CacheA outperforms Sgreedy and Agreedy in terms of total time by 36% and 10% respectively. bCompute stores the cache data on the site with most idle CPUs, which is often the site with the most CPUs. This leads the cached data to be stored close to where it is generated, thus reducing data transfers when adding data to the cache. For the second user, CacheA outperforms Sgreedy and Agreedy in terms of total time by 46% and 21% respectively. The cached data generated by the first user is stored on the sites with more available CPUs, which minimizes the intermediate and reused cached data transfers. From the third user, the storage at some site gets full, i.e. for the third user's execution, Site 3 storage is full and from the fourth user's execution, Site 2 storage is full. Thus, the performance of the three scheduling methods decreases due to higher cache data transfer times. Yet, CacheA still outperforms Sgreedy and Agreedy in terms of total time by 49% and 25% respectively. Experiment 3: cloud site heterogeneity. We now compare the three methods in the case of heterogeneous sites by considering the amount of data transferred and execution time. In this experiment (see Figure 3), we consider only one user who executes the workflow and that previous executions with 60% of the same raw data have generated some cached data. We use the bStorage method for cache site selection. With heterogeneous sites (H > 0), the sites with more CPUs have less available storage but can execute more tasks, which leads to a larger amount of intermediate and cached data transferred between the sites. For H = 0.3, CacheA outperforms Sgreedy and Agreedy in terms of total time (by 40% and 18%, respectively) and amount of data transferred (by 47% and 21%, respectively).

With H = 0.7, CacheA outperforms Sgreedy and Agreedy in terms of total time (by 58% and 42%, respectively) and in terms of amount of data transferred (by 55% and 31%, respectively). CacheA is faster because its scheduling leads to a smaller amount of cached data transferred when reused (48% smaller than Agreedy) and added to the cache (62% smaller than Agreedy).

Concluding Remarks

Our cache-aware scheduling method CacheA always outperforms the two baseline methods (which also benefit from our cache/data selection method), both in the case of multiple users and heterogeneous sites.

The first experiment (with caching) shows that storing and reusing cached data becomes beneficial when 20% or more of the input data is reused. The second experiment (multiple users) shows that CacheA outperforms Sgreedy and Agreedy in terms of total time by up to 61% and 43%, respectively. It also shows that, with increasing numbers of users, the performance of the three scheduling methods decreases due to higher cache data transfer times. The third experiment (heterogeneous sites) shows that CacheA adapts well to site heterogeneity, minimizing the amount of cached data transferred and thus reducing total time. It outperforms Sgreedy and Agreedy in terms of total time by up to 58% and 42% respectively.

Both cache site selection methods bCompute and bStorage have their own advantages. bCompute outperforms bStorage in terms of data transfer time by 13% for the first user and up to 17% for the second user. However, it does not scale with the number of users, and the limited storage capacities of Site 2 and 3 lead to a bottleneck. On the other hand, bStorage balances the cached data among sites and prevents the bottleneck when accessing the cached data, thus reducing re-execution times. In summary, bCompute is best suited for computeintensive workflows that generate smaller intermediate datasets while bStorage is best suited for data-intensive workflows where executions can be performed at the site where the data is stored.

Conclusion

In this paper, we proposed a solution for distributed caching of scientific workflows in a cloud with heterogeneous sites (including on premise servers and shared-nothing clusters). Based on a distributed and parallel architecture, we proposed algorithms for adaptive caching, cache site selection and dynamic workflow scheduling. We implemented our solution in OpenAlea, together with a multisite scheduling algorithm. Using a real data-intensive application in plant phenotyping (Phenomenal), our extensive experimental evaluation using a cloud with three heterogeneous sites shows that our solution can yield major performance gains. In particular, it reduces much execution times and data transfers, compared to two baseline scheduling methods (which also use our cache/data selection method).

Fig. 1 :

 1 Fig. 1: Multisite Workflow System Architecture

  (a) bStorage method (b) bCompute method

Fig. 2 :

 2 Fig. 2: Total times for multiple users (60% of same raw data per user) for three scheduling methods (Sgreedy, Agreedy and CacheA).

Fig. 3 :

 3 Fig. 3: Execution times and amounts of data transferred for one user (60% of same raw data used), on heterogeneous sites with three scheduling methods (Sgreedy, Agreedy and CacheA).

Figure 3

 3 Figure 3 shows the execution times and the amount of data transferred using the three scheduling methods in case of heterogeneous sites. With homogeneous sites (H = 0), the three methods have almost the same execution time. CacheA outperforms Sgreedy in terms of amount of intermediate data transferred and total time by 44% and 26%, respectively. CacheA has execution time similar to Agreedy (3.1% longer). The cached data is balanced as the three sites have same storage capacities. Thus, the total times of CacheA and Agreedy are almost the same.With heterogeneous sites (H > 0), the sites with more CPUs have less available storage but can execute more tasks, which leads to a larger amount of intermediate and cached data transferred between the sites. For H = 0.3, CacheA outperforms Sgreedy and Agreedy in terms of total time (by 40% and 18%, respectively) and amount of data transferred (by 47% and 21%, respectively).With H = 0.7, CacheA outperforms Sgreedy and Agreedy in terms of total time (by 58% and 42%, respectively) and in terms of amount of data transferred (by 55% and 31%, respectively). CacheA is faster because its scheduling leads to a smaller amount of cached data transferred when reused (48% smaller than Agreedy) and added to the cache (62% smaller than Agreedy).

https://www.phenome-emphasis.fr/phenome_eng/

https://cassandra.apache.org

Acknowledgments

This work was supported by the #DigitAg French initiative (www.hdigitag.fr), the SciDISC and HPDaSc Inria associated teams with Brazil, the Phenome-Emphasis project (ANR-11-INBS-0012) and IFB (ANR-11-INBS-0013) from the Agence Nationale de la Recherche and the France Grille Scientific Interest Group.


Sarah Laroui

Xavier Descombes

Aurélia Vernay

Florent Villiers

François Villalba

Eric Debreuve

How to define a rejection class based on model learning?

In supervised classification, the learning process typically trains a classifier to optimize the accuracy of classifying data into the classes that appear in the learning set, and only them. While this framework fits many use cases, there are situations where the learning process is knowingly performed using a learning set that only represents the data that have been observed so far among a virtually unconstrained variety of possible samples. It is then crucial to define a classifier which has the ability to reject a sample, i.e., to classify it into a rejection class that has not been yet defined. Although obvious solutions can add this ability a posteriori to a classifier that has been learned classically, a better approach seems to directly account for this requirement in the classifier design.

In this paper, we propose an innovative learning strategy for supervised classification that is able, by design, to reject a sample as not belonging to any of the known classes. For that, we rely on modeling each class as the combination of a probability density function (PDF) and a threshold that is computed with respect to the other classes. Several alternatives are proposed and compared in this framework. A comparison with straightforward approaches is also provided.

Introduction

In supervised classification, the learning process typically trains a classifier to optimize the accuracy of classifying data into the classes that appear in the learning set, and only them [START_REF] Siegismund | Developing deep learning applications for life science and pharma industry[END_REF]. This framework fits many use cases. For example, in a medical context, the expert knows in advance that the data will be distributed among the healthy class and a given number of pathological classes. There are situations though where the learning process is knowingly performed using a learning set that only represents the data that have been observed so far among a virtually unconstrained variety of possible samples. For example, in a biological study where cells are cultivated in vitro in the presence of a chemical compound, the activity of such molecule can result in specific and remarkable morphological cellular changes collectively called a phenotypic signature. One might want to build a classifier ensuring the automatic recognition of known phenotypes, ensuring in the future a robust and valid clustering of such cases. In parallel, all new chemical treatment could lead to new morphological changes or phenotypes, which should be automatically assigned to the so-called rejection class. Another example is given by a genetic screen study where the different possible patterns induced by inhibiting a given gene are unknown. Obvious solutions can add this ability to a classifier that has been learned classically. For a k-Nearest neighbor classifier, a threshold on the distance to the neighbors can be used to reject samples that are too far [START_REF] Dubuisson | A statistical decision rule with incomplete knowledge about the classes[END_REF][START_REF] Hellman | The nearest neighbor classification rule with a reject option[END_REF][START_REF] Dalitz | Reject options and confidence measures for knn classifiers[END_REF]. For a Deep Neural Network, a function of the output layer can be defined so as to represent a confidence of the prediction. Applying a threshold on such a confidence allows to reject samples [START_REF] Chow | On optimum recognition error and reject tradeoff[END_REF][START_REF] Vasconcelos | Investigating feedforward neural networks with respect to the rejection of spurious patterns[END_REF]. Nevertheless, a classifier is a feature space splitter which focuses on optimizing the frontiers between classes without much attention to the space in-between classes or away from them. Therefore, a better approach seems to directly account for the requirement to possibly reject samples in the classifier design.

In [START_REF] Stefano | To reject or not to reject: that is the question-an answer in case of neural classifiers[END_REF], a rejection option is defined on the output of a Convolutional Neural Networks (CNN). The rejection procedure is optimized based on three costs: the (negative) cost of making a good prediction, the cost of making a wrong prediction, and the cost of rejecting a sample, which is lower than when predicting wrongly. The procedure is optimal when it rejects as much bad predictions as possible while keeping as much good predictions as possible. It relies on thresholds which are estimated using a learning set.

Figure 1: Illustration of the dependency of the FMF thresholds on the proximity between the classes using three Gaussians with the same variance. This paper presents a method for supervised classification with a rejection class which accounts for the neighboring relationships between the known classes instead of classically defining a rejection procedure as a collection of independent, per-class decisions. The main steps are the model estimations for each class and the learning of thresholds based on inter-class relationships. The paper is organized as follows: Section 2 description of the proposed method, the decision rules and the proposals for establishing a model and a threshold specific for each class. In Section 3, a case study in the biochemistry domain and how our method is employed to classify in known class and reject samples which not belonging to any of the known classes. Then we will conclude based on the results obtained with the different methods.

Proposed method

The proposed learning strategy relies on class models where a model is a combination of a fuzzy membership function (FMF) (typically a probability density function (PDF)) on the data space and a threshold. For a sample s, if F is an FMF of a class, then F (s) is a real value between zero and one indicating how likely it is that s belongs to the considered class. The FMFs are learned independently for each class [START_REF] Brew | An evaluation of one-class classification techniques for speaker verification[END_REF]. However, the thresholds are learned by taking into account the overlapping between the classes, thus inferring the possible location of new classes depending on the variability of neighboring classes (see Fig. 1). Finally, the prediction of a sample class is made by considering the responses from all the models. It can either fall into one of the classes used for learning or into the rejection class.

Notations and sketch of the proposed strategy

Let S be the learning set

S = ∪ nc i=1 S i (1)
where S i is the subset of samples of class C i and n c is the number of classes defined by the expert of the domain. The samples of S i are denoted by s i,j where j is the sample index between 1 and n Si . The sketch of the learning procedure is

• For all i, learn the FMF M i of C i based on S i ;

• For all i and for a given k, learn a threshold T i,k on M i based on M i , S i , and ∪ j∈J k S j where J k is defined in Eq. [START_REF] Dubuisson | A statistical decision rule with incomplete knowledge about the classes[END_REF].

J k ⊂ N c = {1, 2, . . . , n c }, i ∈ J k . (2)
Let us consider two options for J k :

-

If J k is N c \ {i}, then the threshold T i on M i is set to T i,k ;
-Otherwise the threshold T i is set to some function of the T i,k 's.

For some sample s, the sketch of the prediction procedure is

• For all i, compute the hard membership response M i (s) ≥ T i ;

• If all responses are false, then s is classified into the rejection class;

• If only the j th response is true, then s is classified into class C j ;

• If several responses are true, then classify s into the class j corresponding to the highest FMF value among the true responses (see Eq. ( 3)).

j = arg max k|M k (s)≥T k M k (s) (3)
Section 2.2 details the proposed model learning approach. Note that Section 2.2.4 describes an alternative to the above prediction procedure.

Model learning

Fuzzy membership function (FMF)

A typical way to learn an FMF for the class C i based on S i is to fit a Gaussian Mixture Model (GMM) [START_REF] Biernacki | Assessing a Mixture Model for Clustering with the Integrated Classification Likelihood[END_REF] to the samples of S i . The number of components of the mixture is a parameter whose optimization is discussed in Section 3.6. For simplicity, it will be chosen identical for all the classes. Moreover, we decided to keep the Gaussian components isotropic. The position and variance of the components will be optimized classically using the expectation-maximization (EM) algorithm [START_REF] Biernacki | Assessing a Mixture Model for Clustering with the Integrated Classification Likelihood[END_REF].

Thresholds on the FMFs

We propose to learn the threshold T i on the FMF M i of the class C i by analyzing how this FMF responds to "its own" sample set S i , naturally, but also to the sample sets of the other classes. There are then two main options. (i) In the 1-vs-all approach, T i is directly learned from S i and all the S k 's, k = i. (ii) In the 1-vs-some approach, a set of thresholds T i,k 's is first learned from S i and some S k 's, k ∈ J k (see Eq. ( 2)). Then the T i,k 's are combined to produce the threshold T i (see Section 2.2.3). One particular case of the 1-vs-some approach is 1-vs-1 where n c -1 thresholds are learned by considering all pairs (S i , S k ), k = i.

For any of the 1-vs-options, the problem is then to learn a threshold T based on M i , S i , and some S k 's, k ∈ J k . For each sample s in the given sample sets, there is a value M i (s). Ideally, this value is high if s is in S i , and low otherwise. Besides, we can symbolically tag the samples of S i with the value one and the samples of the other sets with the value zero. Then, we can plot the tag value against the value of M i for each sample. Once presented this way, an idea that comes to mind to learn T is to perform a logistic regression [START_REF] Sperandei | Understanding logistic regression analysis[END_REF]. The obtained regression sigmoid R moves from the zero-tag position up to the one-tag one. It seems appropriate to define the threshold T on M i as the inverse image of the fictitious tag value 0.5

T = R -1 (0.5). (4)
A potential drawback of this estimation is its sensitivity to outliers. Indeed, depending on the regularization applied to the regression problem, the transition of R from zero to one might be shifted in the presence of outliers (samples of S i with a very low value of M i and samples of some S k , k = i, with a very high value of M i ). The severity of such a shift will depend on the horizontal position of the outliers (see Fig. 2). Thus, we propose an alternative method to learn T using the number of misclassifications. If a sample s belongs to S i , T must be such that M i (s) ≥ T is true. If not, the prediction of s is a false negative. Conversely, the condition M i (s) < T must be true if s belongs to the other samples sets. If not, the prediction is a false positive. The threshold T can then be expressed as the/a minimizer of

G(T ) = #{s ∈ S i s.t. M i (s) < T } -#{s ∈ S k , k ∈ J k s.t. M i (s) ≥ T } (5)
where #S is the cardinal of S. In practice, G does not have a unique minimizer so that T cannot be defined as arg min T G(T ). Instead, there is usually an interval [T min , T max ] such that

∀T ∈ [T min , T max ], G(T ) = min U G(U ). ( 6
)
Figure 2: FMF threshold learning illustrated on synthetic data with the two options: regression-based (top) and misclassification-based (bottom). Each color of curve and dashed thresholding line corresponds to a learning result with the learning set composed of the black cross and dot samples plus the cross sample of the same color. It can be observed that shifting the outlier position also shifts the threshold in the same direction when learning with regression, whereas the threshold learned using the misclassification criterion does not change (as long as the outlier does not enter the [T min , T max ] as mentioned in the text). Note that outliers were assigned a weight of five in order to emphasize their effect. Then, we propose to define T as a function of T min and T max (see Fig. 2). The typical choice is the average function

T = (T min + T max )/2. ( 7
)
Regarding the sensitivity of the estimation (7) of T to outliers, note that as long as no outlier moves into the interval [T min , T max ] (i.e., the values of M i for the outliers remains outside this interval), the estimated value is unchanged (see Fig. 2). However, the presence of an outlier influences the threshold learning in a discrete way, whereas its influence is continuous in the regression-based approach. By discrete, we mean that the presence of an outlier can make T jump from one interval between samples (a [T min , T max ] interval) to a neighboring interval (see vertical dashed lines in Fig. 2). By continuous, we mean that the presence of an outlier will change T in a smoother way.

The pros and cons of the regression-based estimation of T and its misclassification-based estimation are summarized in Table 1.

Combining the FMF thresholds in 1-vs-1 context

As explained in Section 2.2.2, in the 1-vs-1 context, n c -1 thresholds T i,k are learned by considering all pairs (S i , S k ), k = i. The threshold T i must then be computed as a function of the T i,k 's. Three straightforward options are the minimum, average, and maximum. Selecting the minimum among the T i,k 's means that T i is chosen according the class that is the farthest to C i . This would increase prediction errors for samples of the closer classes. Therefore, we propose to defined T i as the maximum of the T i,k 's, thus limiting the confusion between nearby classes.

Classification into known classes

The classification procedure into the rejection class has been defined in Section 2.1. The model-based procedure for classifying into known classes has also been described in the same section. We now suggest an alternative using an external classifier. Indeed, since classifiers such as CNN often exhibit very good classification performances in the no-rejection-class framework, they might be used advantageously in replacement of the model-based procedure whenever a sample has not been rejected (at least one hard membership response is true). Using a CNN would also provide us with the feature representation for our data.

3 Experimental results

Biological context

Botrytis cinerea is a reference model of filamentous phytopathogen fungi. We have observed that some anti-fungal molecules associated to an hypothesis of Mode of Action (MoA) can lead to characteristic morphological changes, or phenotypic signatures (see Fig. 3). In order to recognize and characterize the different phenotypes and associate them with the various modes of action, images are acquired by transmitted light microscopy after 24 hours of culture in the presence of molecules. Then, the images are labeled by an expert according to the observed phenotype.

Feature vectors for model estimation

Deep Neural Network (DNN)

Instead of hand-crafting features, we decided to take advantage of the feature learning capability of Convolutional Neural Networks (CNN). As mentioned in Section 2.2.4, we will also used the CNN for its predictions on non-rejected samples. For the sake of simplicity, we chose to rely on the transfer learning [START_REF] Shin | Deep convolutional neural networks for computer-aided detection: Cnn architectures, dataset characteristics and transfer learning[END_REF] property of DNNs by using MobileNet [START_REF] Howard | Mobilenets: Efficient convolutional neural networks for mobile vision applications[END_REF] pre-trained on the ImageNet database [START_REF] Deng | ImageNet: A Large-Scale Hierarchical Image Database[END_REF]. The feature learning layers were frozen while the classification part was removed, replaced with one corresponding to our classification problem into four classes, and trained on our data. The raw data were 2160x2160-pixel grayscale images, collectively denoted by D. Each image was divided into 16 sub-images of 500x500 pixels, then rescaled to 224x224 pixels to fit the input shape of MobileNet. The rescaled sub-images of an image I are denoted by

I k , k ∈ [1..16].

Selected features

Once MobileNet re-trained on our data (see Section 3.2.1), each sub-image could be represented by the feature vector produced by MobileNet from its rescaled version. Let us see now the two proposed options for feature vector definitions.

Bottlenecks

The term "Bottleneck" classically refers to the output of the last layer before the classification part of the network usually made of one or more dense layers. This is the straightforward option for feature vector definition. In MobileNet, this output is a 1001-dimensional vector, which we believe is too high to properly learn our class models. Therefore, a dimension reduction was performed.

Two approaches were evaluated: Principal Component Analysis (PCA) [START_REF] Jolliffe | Principal Component Analysis, ser. Springer Series in Statistics[END_REF] and an approach based on variable importance (VI) as defined in supervised classification [START_REF] Criminisi | Decision forests for classification, regression, density estimation, manifold learning and semi-supervised learning[END_REF]. The latter one consists in computing the variable importance measure of each component of the feature vector definition. This is conveniently provided by publicly available Random Forest libraries [START_REF] Pedregosa | Scikit-learn: Machine learning in Python[END_REF], for example. The features were ranked from highest to lowest importance, and the k most important features were retained for further processing. Classically, the same procedure was applied to PCA-transformed features where variance plays the role of variable importance. In both case, the number k of retained features is a parameter that must be optimized.

Output layer Here, we suggest that the vector of the output layer of a DNN can be viewed as a feature vector instead of viewing it as an informative classification response from which to compute the prediction. Each component reflects a form of confidence that a sample belongs to the corresponding class. It has a value between zero and one, and the sum of all the components is normalized to one. Whenever one component is very close to one (while the others are therefore close to zero), the classification response has only a single reasonable interpretation. However, when there is no clear maximum among the components, the way the component "budget" of one is distributed among them might be interpretable by another classifier in order to correctly classify the corresponding sample. The use of the output layer as a feature vector follows from this hypothesis.

Let us study the vectors of the output layer for all the samples of a class. The histogram of the values of the component associated with the correct class will certainly have high values toward one and low ones toward zero. The opposite happens for the components associated with the other classes. In either case, the histogram looks like a truncated version of a hypothetical symmetric histogram centered on one or zero. As a consequence, if these vectors were to be used as-is as samples to learn a GMM, the mismatch between the non-symmetric component distribution and the symmetry of the Gaussian basis functions would certainly degrade the quality of the model fitting (see Fig. 4). In order to improve the fitting quality, we propose to symmetrize the component distribution by sample duplication. Let s be a sample s = (s 1 , s 2 , . . . , s d ) T .

If the first component s 1 is greater or equal to 0.5, then a new sample with a component s i = 2 -s i is created. If s i < 0.5, then the component of the new sample is equal to -s i . This procedure is repeated for both the sample and the newly created one with the second component. And again until the component s d . At each step, the number of samples is multiplied by two. In the end, the symmetrized sample set contains 2 d more samples then originally. In our case, the output vectors are 4-dimensional, so the symmetrized set has 16 more samples. Note that once the model fitting is done (see Fig. 4), the samples created for symmetrization are simply discarded and the learned model can be used normally on the original samples.

Learning set for model estimation

In practice, we noted that, because of the sparse and non-uniform spatial distribution of the fungi, some sub-images were empty while others were densely populated with objects. The feature vectors of these sub-images correspond to outliers in the learning set. Therefore, we moved from a learning set at the level of sub-images to an image-level learning set

L = {m(f (I 1 ), f (I 2 ), . . . , f (I 16 )) | I ∈ D} (9)
where f denotes the mapping between a rescaled sub-image and its feature vector as produced by Mo-bileNet, and m is the vector of per-component medians of its (vector) arguments (normalized between 0 and 1). The number of samples per known class in L was 233, 225, 297 and 441 for classes 1 to 4 respectively. Classically, a part of L will be used for training, a part will be used to assess the accuracy of the proposed method on the known classes in a cross-validation framework. Moreover, 638 samples that were identified by an expert as not pertaining to any of the four known classes will be used to assess the ability of the proposed method to classify into the rejection class.

Model estimation

As mentioned in Section 2.2.1, the learning of the FMFs of the class models was performed using a GMM fitting. A Gaussian is characterized by its mean and covariance matrix. Three constraint levels on the covariance matrix were tested: spherical (the covariance matrix is diagonal-constant), diagonal (the covariance matrix is diagonal), and full (no constraint on the covariance matrix). Each component of the mixture has its own covariance matrix. If k is the number of components selected for the fitting and d is the number of dimensions of the sample space, then the number of GMM parameters to estimate is k(d + 1) for spherical, k(d + d) for diagonal, and k(d + d(d + 1)/2) for full.

CNN-based rejection class

Method CNNro

The proposed method was compared to the method in [START_REF] Stefano | To reject or not to reject: that is the question-an answer in case of neural classifiers[END_REF], denoted here by CNNro. This work defines a rejection option that can be applied to any classical classifier. The rejection option is based on a function defined to estimate the classification reliability of a sample. Then a threshold is learned on this function so as to reject as much wrongly classified samples as possible while keeping as much correctly classified samples as possible. This is expressed as the optimization of a performance measure P with respect to a threshold σ. In practice, the authors propose to implement their approach for a CNN and use two sets of reliability function and threshold in cascade. The first reliability function, Ψ A , corresponds to the maximum of the output layer. The second function, Ψ B , corresponds to the difference between the two largest values of the output layer. A sample s is rejected if Ψ A (s) < σ A , or else if Ψ B (s) < σ B , where σ A and σ B are optimized according to P on the whole sample set and on the samples that have not been rejected by Ψ A , respectively. This optimization depends on a parameter called the normalized cost C N (see Section 3.6).

Methods CNNa and CNNr

We also compared the proposed method with two straightforward ways to reject a sample with a CNN. (i) If the highest value among the output neurons is lower than some threshold, then the sample is classified into the rejection class. Otherwise, it is classified into the class corresponding to the neuron as usual. Let us call this procedure the absolute-threshold rejection method and denote it by CNNa. Alternatively, (ii) the thresholding can be made on the difference between the two highest values among the output neurons. This procedure is called relative thresholding and denoted by CNNr. In both cases, a threshold T CN N must be learned, being either "a" or "r". In order to perform this learning in the same conditions as the ones for the proposed method, the samples of the unknown class were not used. Instead, each one of the four known classes was used in turn as if unknown (say class C 4 ), the three others playing their actual role of known class (C 1 , C 2 and C 3 ). This means that a new CNN classifier has to be learned for just the current three known classes. It was done in a transfer learning mode: keeping the feature learning layers, cutting out the classification layers, adding new classification layers with three output neurons, and optimizing only this newly added part of the CNN. Afterwards, a threshold T CN N (1, 2, 3|4) was estimated in a brute force way so as to maximize the weighted average accuracy over the known and rejection classes (weight: relative number of samples per class as usual). In practice, a hundred thresholds between 0.3 and 1 were tested. In the end, four thresholds were learned: T CN N ({1, 2, 3, 4} \ {i}|i) where i ∈ [1..4]. The threshold T CN N was computed as a function of the

T CN N (• • • |•).
In our application, the median led to the best results (where the median of four values is simply the average of the two middle ones).

Results and analysis

The experimental results were averaged over a repeated random sub-sampling cross-validation of 20 rounds. In each round, the data corresponding to known classes were randomly split into a training set (50%), a validation set used for parameter optimization (30%), and a test set for the estimation of classification accuracies (20%). In each training, validation and test set, we have respected the balance between the classes. Note that the data corresponding to unknown classes were used only for testing. As a reminder, the images were represented by the bottleneck features or the output of a CNN. When dealing with bottleneck features, dimension reduction was performed using variable importance or PCA, and the optimal reduced dimension was selected using brute force search in the interval [START_REF] Sperandei | Understanding logistic regression analysis[END_REF]300] with a step of 10 from 10 to 100 and a step of 100 afterwards. When dealing with the CNN outputs, the data were symmetrized or not. The FMFs of the classes were estimated using GMMs. We tested different types of covariance parameters, we have kept the ones that give the best results. Note that the number of Gaussian components was optimized globally for all the classes. The optimal number was selected using brute force search in the interval [START_REF] Siegismund | Developing deep learning applications for life science and pharma industry[END_REF][START_REF] Shin | Deep convolutional neural networks for computer-aided detection: Cnn architectures, dataset characteristics and transfer learning[END_REF]. The thresholds on the FMFs were estimated using either logistic regression or a misclassification score, in a 1-vs-all or 1-vs-1 context. The classification accuracies obtained in all these experimental settings are presented in the table 2.

Methods

Features Results

Threshold learning Regression

Misclassif. 1-vs-all 1-vs-1 1-vs-all 1-vs-1 Table 2: Classification accuracies as percentages for the various experimental settings. The actual values were rounded to the closest integer. In particular, values of 100 in the table correspond to accuracies above 99.5%. The "Known" rows correspond to the accuracy on the known classes and the "Rejection" rows correspond to the other, unknown ones (the rejection class from the classifier point of view). A parameter of the GMM corresponds to the type of covariance used. Diagonal is the case where the contour axes are oriented along the coordinate axes and the eccentricities may vary between components. Spherical is a "diagonal" situation with circular contours. CNNro corresponds to the method proposed in [START_REF] Stefano | To reject or not to reject: that is the question-an answer in case of neural classifiers[END_REF] (see Section 3.5.1). CNNa denotes the absolute-threshold rejection method for a CNN. CNNr denotes the method of rejection based on relative thresholding. CNNmax refers to MobileNet learned classically on the four known classes. Therefore, it has no sample rejection capability.

They were computed separately for the test samples of the known classes (see "Known" rows) and the samples of the other, unknown classes (which, as mentioned earlier, were never used in training or validation; See "Rejection" rows). As previously explained, the classification into known classes can be done by the proposed method or falling back to an additional classifier, here the CNN used to compute the feature vectors of the data (see Sections 2.2.4 and 3.2.1).

When using the bottleneck features with the proposed method, the best results are always obtained when reducing the feature space dimension based on variable importance. Within this context, the best results are always obtained with the regression-based threshold learning. Within this "sub-context", the best results are obtained in the 1-vs-1 context. When using the output features, the best results are always obtained when symmetrizing the data (which allows a better GMM fitting). Within this context, the best results are always obtained with the regression-based threshold learning. Within this "sub-context", the best results are obtained in the 1-vs-all context. Overall, the best setting corresponds to the 1-vs-all regression-based threshold learning on FMFs estimated on symmetrized output features. The accuracies were 93% and 92% on the known samples and the samples of the new classes, respectively.

The proposed method was compared with the method CNNro (see Section 3.5.1). To find the two thresholds σ A and σ B , we tested a range of values for the normalized cost C N , which is given by

C N = C e C r -1 (10)
where C e is the cost of a single misclassification and C r is the cost of a rejection. Using cross-validation on 50 folds, we obtained the best classification accuracies on known and new class samples for C N equal to 25. The optimal values for σ A and σ B were 0.94 and 0.91, respectively. These accuracies were 65% on the known samples and 100% on the samples of the new classes (values rounded to the closest integer). The proposed method is largely better on the known samples (+28%). It performs worse on the new samples, but not dramatically (92% versus 100%). Together with the fact that the proposed method is also more balanced (1% difference between known and new versus 35% for CNNro), it can be concluded that it is significantly better. The proposed method was also compared with the CNN with absolute threshold method (CNNa) and the CNN with relative threshold method (CNNr) (see Section 3.5.2). In the learning stage, the following thresholds were found: {0.70, 0.72, 0.84, 0.98} for T CN N a (• • • |•), and {0.49, 0.52, 0.73, 0.96} for T CN N r (• • • |•). In both cases, there is a significant variability among the thresholds. It probably reveals that the four known classes do not overlap each other equally. A tentative explanation goes as follows: if the rejection class is far away from the three known classes, the CNN output for its samples will certainly have fairly low values for all of the output neurons. The rejection procedure can then afford a low threshold. On the contrary, if the rejection class overlaps significantly the three known classes, the CNN output for its samples will certainly be closer to the outputs for the samples of the known classes, that is, one output neuron with a value close to one, the others with a value close to zero. In this case, the rejection procedure is forced to use a high threshold. In practice, taking the median of these thresholds led to the best results. So we used T CN N a = 0.81 and T CN N r = 0.675. The proposed method with its best settings compares very favorably to these two straightforward approaches to defining a rejection class with a CNN classifier. Indeed, it provides accuracies higher by 10 to 25% for both the rejection and the known classes. Besides, it is even on par with the (unmodified) learned CNN (denoted by CNNmax) on the known classes with an accuracy of 93% versus 91% for CNNmax (see Table 2).

Conclusion and perspectives

We proposed a method of rejection class based on model learning in a supervised context. It follows a general strategy with three main steps: model learning independently for each class, threshold learning based on class interactions, and prediction procedure. Although the focus is on sample rejection, the method can also be used for classification into the classes known at the time of learning.

In the studied application, the proposed method provided very good classification accuracies for both the rejection class and the known classes. On the known classes, it proved to be equally effective than a CNN classifier. The significantly lower performances of two straightforward approaches to defining a rejection class with a CNN classifier tend to show that defining a rejection class should probably not be considered as an afterthought, but rather be part of the classifier design.

Regardless of the performances, let us note the main practical differences between bottleneck and output features. With bottleneck, reducing the feature space dimension will often be necessary. However, the features do not have to be learned again if a new class is formed. With the output layer, the feature space dimension is reasonable. Nevertheless, adapting and re-training the classification part of the CNN is mandatory whenever new classes enrich the expert knowledge.

As noted above, the performances in our application are very satisfying. However, the accuracy is slightly better for samples of known classes with the bottleneck features. Conversely, the accuracy is better for rejected samples with output features (except the notable "Original/Regression/1-vs-all" case). Future works will concentrate on understanding and reducing or removing these biases.

Experimentally, there is an improvement we would like to consider. Here, the model fitting to learn the FMFs was done in the same conditions for all the known classes. Namely, the GMMs were optimized for the same number of components for all the classes. Letting the number of components be also optimized on a per class basis could improve the results further. In particular, the number of components could be estimated using an RJMCMC algorithm [START_REF] Green | Reversible jump Markov chain Monte Carlo computation and Bayesian model determination[END_REF], although it would imply a high computational time.

Figure 3 :

 3 Figure 3: Four different classes observed on Botrytis cinerea model: four known characteristic phenotypes associated with different modes of action linked to different antifungal chemical families.

Figure 4 :

 4 Figure 4: Fitting a Gaussian onto the histograms of the four features for the samples of one given class: (top row) original data and (bottom row) symmetrized data. The Gaussians better approximate the original data (abscissa between 0 and 1) when fitted onto the symmetrized data.

Table 1 :

 1 Sensitivity to outliers of the regression-based and the misclassification-based estimations of T (see text for details).

			=Outliers
		Influenced by position	Influence per
	Regression	Yes	Continuous
	Misclassif.	No	Discrete