Leo Wanner
email: leo.wanner@upf.edu

Stefanos Vrochidis
email: stefanos@iti.gr

Marco Rospocher
email: rospocher@fbk.eu

Jürgen Mossgraber
email: juergen.mossgraber@iosb.fraunhofer.de

Harald Bosch
email: harald.bosch@vis.uni-stuttgart.de

Ari Karppinen
email: ari.karppinen@fmi.fi

Maria Myllynen
email: maria.myllynen@hsy.fi

Sara Tonelli
email: satonelli@fbk.eu

Nadjet Bouayad-Agha
email: nadjet.bouayad@upf.edu

Gerard Casamayor
email: gerard.casamayor@upf.edu

Jürgen Moßgraber

Thomas Ertl
email: thomas.ertl@vis.uni-stuttgart.de

Désirée Hilbring
email: desiree.hilbring@iosb.fraunhofer.de

Lasse Johansson
email: lasse.johansson@fmi.fi

Kostas Karatzas

Ioannis Kompatsiaris

Tarja Koskentalo
email: tarja.koskentalo@hsy.fi

Simon Mille
email: simon.mille@upf.edu

Anastasia Moumtzidou
email: moumtzid@iti.gr

Emanuele Pianta
email: pianta@fbk.eu

Luciano Serafini
email: serafini@fbk.eu

Virpi Tarvainen
email: virpi.tarvainen@fmi.fi

Kessler

Personalized Environmental Service Orchestration for Quality of Life Improvement

Keywords: environmental information service, environmental node discovery, knowledge, personalization, infrastructure, services

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

Introduction

Environmental and meteorological conditions are considered of utmost importance for the population, as they strongly influence the quality of human life. Citizens are increasingly aware of the important role that environmental data and measurements play in health issues (e.g. allergies), as well as in a variety of important daily activities (e.g. agriculture). One of the consequences of this awareness is the demand for high quality environmental information that is tailored to one's specific context and background (e.g. health conditions, travel preferences, etc.), i.e., which is personalized. Personalized environmental information may need to cover a variety of aspects (such as meteorology, air quality and pollen) and take into account a number of specific personal attributes (health, age, etc.) of the user, as well as the intended use of the information. So far, only a few approaches have been proposed with a view of how this information can be facilitated in technical terms. All of these approaches focus on one environmental aspect and only very few of them address the problem of information personalization [START_REF] Karatzas | State-of-the-art in the dissemination of AQ information to the general public[END_REF], [START_REF] Peinel | Customized Information Services for Environmental Awareness in Urban Areas[END_REF], [START_REF] Wanner | MARQUIS: Generation of User-Tailored Multilingual Air Quality Bulletins[END_REF]. On the contrary, we aim at addressing the aforementioned task in its full complexity.

In our work, we take advantage of the fact that nowadays, the World Wide Web already hosts a great range of services (i.e. websites, which provide environmental information) that offer data on each of the above aspects, such that, in principle, the required basic data are available. The challenge is threefold: first, to discover and orchestrate these services; second, to process the obtained data in accordance with the needs of the user; and, third, to communicate the gained information in the user's preferred mode. To address this problem, we need to involve a considerable number of rather heterogeneous applications and thus create an infrastructure that is flexible and stable enough to support a potentially distributed architecture. This infrastructure is realized by the "PESCaDO platform", which is described in the upcoming sections.

In Section 2, we outline the process of the discovery and extraction of environmental information from services (also referred to as nodes) in the Web, which is considered as the prerequisite step for enabling the retrieval capabilities of the system. In Section 3, we describe briefly the processing of the data obtained from the environmental nodes until their delivery to the user. In Section 4, we present the infrastructure designed to accommodate for both the discovery itself and the posterior tasks. Section 5 illustrates the functionality of our system, while Section 6 concludes with a brief summary of our exposé.

A Sample User Scenario

In order to set up the context of the use of the PESCaDO platform, let us imagine a non-professional user, Fiona Fit, who is in her late twenties and lives in the municipality of Espoo, located in the south area of Finland. Fiona is rather active in her leisure time and often goes hiking. But since she is allergic to birch pollen, she needs information on the environmental conditions in the hiking resorts before she decides on her route. This afternoon, Fiona wants to go for a hike in the Nuuksio National Park near Helsinki and needs to know whether the forecasted air quality, weather and pollen conditions are favorable. As a registered user, with her profile uploaded to the system, she seeks decision support from PESCaDO.

Fig. 1. PESCaDO demonstrator and input for the sample use case

Fiona uses PESCaDO's interface to formulate the aforementioned query. 1 Figure 1 displays the interface and the type of input information a user can submit, i.e., the profile, the request type (whether it is an instruction, report or warning request), the type of activity (e.g. travelling, physical outdoor activity), the start/ end date/ time and the region of interest (depicted as blue polygon on the map).

Based on its knowledge regarding the forecasted air quality, pollen and weather values as well Fiona's health status, the system provides an answer which discourages Fiona from engaging in any sport activities because high concentration of Thoracic Particles is expected during the selected time in the region of Nuuksio Park. The system's answer is displayed in Figure 2.

Although a higher degree of personalization still seems possible, especially in terms of a direct amicable salutation of the user and the like, the offered personalized information without any doubt is already capable of improving the quality of life of their addressees.

Let us discuss in the following sections how the system processes a query and how it comes up with a personalized answer. 1 Link to the demonstrator is available at: http://www.pescado-project.eu/

Architecture of the PESCaDO Platform

In order to be able to "understand" the query of the user, to gather the required data, to interpret them and then to generate a recommendation or any other type of information useful to the user, PESCaDO needs to tackle a number of tasks:

1. Discovery of the relevant environmental service nodes in the web and extraction of information. As already mentioned above, the web hosts a large amount of environmental (meteorological, air quality, pollen, etc.) distributed services, which include both public webpages that offer environmental information worldwide, as well as dedicated environmental web services with free access; especially the number of meteorological services that cover each major location is impressive. Among all these services those must be discovered that maybe of relevance to PESCaDO and the information offered by them must be extracted. The heterogeneous forms and formats, including text and images, in which this information is encoded make the task of discovery and extraction of information from webpages that provide environmental information a serious challenge. The information extracted from the discovered service nodes is stored in a repository and indexed (together with the references to the nodes). This task is performed off-line, i.e., independent of the queries of PESCaDO's users, while all the other tasks are performed on-line.

Identification of user relevant service nodes.

The indexed environmental repository compiled as result of the node discovery task and updated at a predefined time rate serves as the general data source for PESCaDO. That is, when a user poses a query, a process of the identification of environmental service nodes in the compiled repository that are relevant to the query of the user, their profile and their context must be launched. This is not trivial, given that, for instance, a user may be moving, be located in an area which is not directly covered by any node.

3. Orchestration of environmental service nodes: Environmental nodes may provide competing or complementary data on the same aspect for the same or the neighboring location. To ensure the availability of the most reliable and most comprehensive content, the contents proceeding from these nodes must be (i) assessed with respect to its trustworthiness and certainty and selected accordingly (if several nodes offer competing data) or (ii) fused (if several nodes offer complementary information.

4. Converting the data into content. In order to guarantee a motivated orchestration of heterogeneous environmental service nodes and offer user-tailored decision support services and environmental information production, we need to convert the data into structured unified content, which allows for application of intelligent reasoning algorithms. To this end, the extracted and fused environmental data are integrated into an environmental knowledge base (KB). Our KB, which is codified in the standard semantic web ontology language OWL [START_REF]World Wide Web Consortium: OWL Web Ontology Language Reference[END_REF], covers environmental content such as meteorological conditions and phenomena, air quality, and pollen, as well as other relevant environment-related content essential for the targeted user-tailored service: travel information, human diseases, geographical data, user profile, etc. In addition, the KB is also capable of formally representing the description of the user's inquiry. The current version of the KB contains around 202 classes, 143 attributes and properties, 463 individuals2 . Its Description Logic expressivity is ALCHOIQ(D). The KB has been obtained by (i) including customized version of currently available ontologies (e.g., parts of the SWEET ontology), (ii) automatically extracting key concepts from domain relevant text sources, and (iii) manually adding additional properties and attributes.

5. Assessment of the content with respect to the needs of the user: Once the data from the nodes have been incorporated into the KB, they need to be evaluated and reasoned about in order to infer how they affect the addressee, given his/her personal health and life circumstances and the purpose of his/her request. For instance, a citizen may request information, because he/she wants to decide upon a planned action, be aware of extreme episodes or monitor the environmental conditions in a location. [START_REF] Boser | A training algorithm for optimal margin classifiers[END_REF]. Selection of user-relevant content and its delivery. Not all content in the KB is apt to be communicated to the addressee: some of it would sound trivial or irrelevant. Intelligent content selection strategies take into account the background of the user and the intended use of the information to decide which elements of the content are worth and meaningful to be communicated. To deliver the selected content, techniques are required that present the content in a suitable mode (text, graphic and/or table) and in the preferred language.

7. Interaction with the user. The design and development of intuitive interfaces for the interaction between the system and the user forms the final task. The user must be able to formulate the problem in a simple and intuitive format and receive the generated information in a suitable form.

In order to accommodate for all tasks described above, we opted for a service-based architecture. This architecture is based on a methodology which has been developed in ORCHESTRA [START_REF]Reference Model for the ORCHESTRA Architecture Version 2.1[END_REF] for risk management, and which has been extended in SANY [START_REF] Usländer | Specification of the Sensor Service Architecture[END_REF] to cover the domain of sensor networks and standard-based sensor web enablement. The focus of this methodology is on a platform-neutral specification. In other words, it aims to provide the basic concepts and their interrelationships (conceptual models) as abstract specifications. The design is guided by the methodology developed in the ISO/IEC Reference Model for Open Distributed Processing (RM-ODP), which explicitly foresees an engineering step that maps solution types, such as information models, services and interfaces specified in information and service viewpoints, respectively, to distributed system technologies.

We defined application-specific major tasks and actions as abstract service specifications, which can be implemented as service instances on a specific platform. Web service instances for these services are currently being developed. They can be redefined and substituted in the future as needed. Figure 3 displays a simplified sample workflow with the major application services in action. Two services are not cited in Figure 3 since they are consulted by nearly all other services: the Knowledge Base Access Service and the User Profile Management Service. The figure also does not include the services related to data discovery and information distillation from webpages.

A main dispatcher service (called Answer Service, AS) controls the workflow and the execution of the services. First, the user interacts with the system via the User Interaction Service (UIS). In the case that the user is unsure with respect to the types of information they can ask for, they can inquire this information by requesting it from the Problem Description Service (PDS).

To ensure a full comprehension of the problem or user generated question, we decided to operate with controlled graphical and natural language input formats. Once the user has chosen what kind of question they want to submit to the system, the UIS provides them the corresponding formats. Thereupon, the user can formulate their query, which is subsequently translated by the PDS into a formal ontology-based representation understood by the system. After the problem description is generated, it is passed by the UIS to the AS as a "Request Answer" inquiry. Then, the AS assesses what kinds of data beyond environmental data are required to answer the query of the user and solicits these data from the Auxiliary Services (AuxS). Such services can provide for instance travel route information in the case that the user's query concerns the environmental conditions for a bicycle tour from A to B.

After having acquired the complementary data, the AS can request from the Data Retrieval Service (DRS) the environmental data needed to answer the user query. The DRS solicits these data from the environmental nodes identified by the Data Node Retrieval Service (DNRS), which accesses the data node repository.

As already mentioned, the retrieved nodes may deliver complementary or competing data of varying quality (to keep the presentation simple, we dispense with the illustration of the orchestration of service nodes). The Fusion Service (FS) applies uncertainty metrics to obtain the optimal and maximally complete data set, which is passed by the AS to the Decision Service (DS). The DS converts the data set into knowledge, in that it relates it to the knowledge in our KB, reasons about it, and assesses it from the perspective of its relevance to the user. From this content, the Content Selection Service (CSS) compiles a content plan which contains the knowledge to be communicated to the user as the answer. The Information Production Service (IPS) takes the content plan as input and generates information in the language and mode (text, table, or graphic) of the user preference, which then is passed to the user.

Fig. 3. Sequence diagram service execution for delivery of environmental information

A number of the above services as well as the interaction between selected services have already been discussed in other publications; see, for instance, [START_REF] Epitropou | Fusion of Environmental Information for the Delivery of Orchestrated Services for the Atmospheric Environment in the PESCaDO project[END_REF] for the presentation of environmental node orchestration in PESCaDO; [START_REF] Rospocher | Ontology Management in a Service-oriented Architecture[END_REF] for the ontology management and [START_REF] Bouayad-Agha | From Ontology to NL: Generation of Multilingual User-Oriented Environmental Reports[END_REF] for the interaction of the DS and CSS and IPS. Therefore, and also due to the lack of space, let us focus in what follows on one aspect of PESCaDO -namely the discovery and extraction of environmental information from the web.

Discovery and Extraction of Environmental Information

The discovery of environmental nodes can be considered a problem of domainspecific web search and therefore methodologies from this area can be applied to implement a node discovery framework; see Figure 4 for the architecture. We apply two types of methodologies of domain search: (a) the use of general purpose search engines for the submission of domain-specific queries, and (b) focused crawling of predetermined websites [START_REF] Wöber | Domain Specific Search Engines[END_REF]. To generate queries for the general purpose search engine we combine domain information from ontologies and geographical data obtained by geographical web services. In this case, the Yahoo BOSS API 3 is used as a general purpose search engine. The queries are expanded by keyword spices [START_REF] Oyama | Domain-Specific Web Search with Keyword Spices Awareness in Urban Areas[END_REF], which are domain specific keywords extracted with the aid of machine learning techniques from environmental websites. In parallel, a set of predefined environmental websites is enriched using a focused crawler, which is capable of exploring the web in a directed fashion in order to collect other nodes that satisfy specific criteria related to the content of the source pages and the link structure of the web. The focused crawler is built upon Apache Nutch ( http://nutch.apache.org/) Fig. 4. Architecture for the discovery of environmental service nodes Since the output of the search engine and the focused crawler include also many irrelevant results we employ a post-processing classification step in order to improve the precision of the discovery phase, which is realized with the aid of Support Vector Machine (SVM) [START_REF] Boser | A training algorithm for optimal margin classifiers[END_REF] classifiers trained with manually annotated shots and textual features extracted with KX [START_REF] Pianta | KX: A Flexible System for Keyphrase Extraction[END_REF], which is a key phrase extraction tool.

The whole discovery procedure is automatic. However, an administrative user could intervene through an interactive graphical user interface in order to select geographic regions of interest to perform the discovery, optimize the selection of keyword spices and parameterize the training of the classifiers. Figure 3 above shows the architecture of the discovery of environmental service nodes.

The next step is to extract the measurement data from the environmental nodes in order to index them in a database. It was observed that usually the weather information is mostly encoded in textual formats in the websites, while the air quality and pollen are usually presented in heatmap images.

To distill the textual data from webpages, advanced natural language processing techniques are needed for webpage parsing, information extraction and text mining. Although these techniques can be tuned to deal with the presentation mode of environmental data and information, this task remains very challenging and still only a manually assisted extraction of such information can be supported due to the high variety of websites and presentation formats. In the case of image data, also a semiautomatic procedure is realized. Specifically, we have implemented a configuration tool, in which the administrative user can identify the important images after the discovery process and annotate specific areas in the heatmap (i.e. coordinates, title, map, etc.). Once the configuration is finalized, the template file is used by a data extraction service to automatically retrieve information from this image in the future. The extraction service is built upon OCR technologies for identifying text on the image, while the image is converted to numerical data with the AirMerge tool [START_REF] Epitropou | A method for the inverse reconstruction of environmental data applicable at the Chemical Weather portal[END_REF].

Finally the environmental node data are stored and indexed in a Sensor Observation Service (SOS) [START_REF]Sensor Observation Service (SOS)[END_REF] compliant repository.

For further details on the discovery of environmental nodes in PESCaDO, see [START_REF] Moumtzidou | Discovery of Environmental Nodes in the Web[END_REF].

Conclusions

In this paper, we have presented a personalized platform based on heterogeneous technologies which can support individuals in planning activities with respect to environmental conditions. As shown in the demonstration, the system can improve the quality of life, since it can offer very important suggestions to people taking into account their health conditions, the indented activities and the environmental conditions.

Fig. 2 .

 2 Fig. 2. PESCaDO demonstrator and output for the sample use case

These data refer to the "empty" KB, i.e. without considering any environmental data coming from the nodes.

Acknowledgments. This work is partially funded by the European Commission under the contract number FP7-248594 "Personalized Environmental Service Configuration and Delivery Orchestration" (PESCaDO).


Shriprakash Sinha
email: shriprakash.sinha@gmail.com

Online Cluster Approximation Via Inequality

Keywords:

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

Introduction

This work focuses on approximating the number of clusters in an online unsupervised learning framework with no prior knowledge of the number of clusters. The current work achieves this using a transductive-inductive (TI) learning strategy. Motivation for the proposed work comes from the recent works on conformal learning theory (CLT) [START_REF] Vovk | Algorithmic learning in a random world[END_REF] and [START_REF] Shafer | A tutorial on conformal prediction[END_REF], where TI has been widely used to estimate the prediction of an unlabeled example based on the already processed data. Also, the estimation of clusters using a multidimensional data is an open area of research. Finding motivation from [START_REF] Kanungo | An efficient-means clustering algorithm: Analysis and implementation[END_REF] on k-means clustering using the Euclidian ball multivariate Chebyshev inequality (MCI) and [START_REF] Sinha | Bounded multivariate surfaces on monovariate internal functions[END_REF] on image representation via a hybrid model based on space filling curve, the proposed work exploits a generalization of MCI [START_REF] Berge | A note on a form of tchebycheff's theorem for two variables[END_REF], to capture the interactions among examples with multiple features. MCI provides bounds for multidimensional data which is afflicted by the curse of dimensionality that make estimation of multivariate probabilities a difficult task. A generalization of MCI is the consideration of probability content of a multivariate normal random vector to lie in an Euclidean n-dimensional ball [START_REF] Monhor | Understanding the concept of outlier and its relevance to the assessment of data quality: Probabilistic background theory[END_REF]. This work employs a conservative approach of using the Euclidian n-dimensional ellipsoid to restrict the spread of the probability content [START_REF] Chen | A new generalization of chebyshev inequality for random vectors[END_REF]. In abstract terms, let X be a stochastic variable in N dimensions with a mean E[X]. Further, let Σ be the covariance matrix of all examples, each containing N features and C p ∈ R (MCI parameter), then the MCI in [START_REF] Chen | A new generalization of chebyshev inequality for random vectors[END_REF] states that:

P{(X -E[X]) T Σ -1 (X -E[X]) ≥ C p } ≤ N
Cp and can be transformed into

P{(X -E[X]) T Σ -1 (X -E[X]) < C p } ≥ 1 -N Cp i.e
. the probability of the spread of the value of X around the sample mean E[X] being greater than C p , is less than N /C p .

The novelty of the current work is in combining the TI learning strategy with the MCI, as the algorithm processes the examples in an online unsupervised paradigm. Before presenting the details, a succinct view of why and how the task is accomplished is given: • The algorithm aims to reconstruct the information context by building clusters, as it processes a random sequence of unlabeled examples. The number of clusters is not known apriori. • A compromise level is used to decide how much loss of information can be tolerated by imposing constraint on the homogeneity of the generated clusters. The compromise level is captured by C p . C p shares similar but not the same concept of defining confidence level as ε does in the CLT paradigm. The level of confidence generated in the proposed algorithm is 1 -N Cp as opposed to 1 -ε in CLT. • The quality of the prediction in CLT is checked based on the p-values generated online. The proposed algorithm generates reconstruction error online as a measure of the quality of reconstruction. The reconstruction error is computed as follows: [START_REF] Vovk | Algorithmic learning in a random world[END_REF] All processed examples present in a cluster are assigned mean feature values by averaging examples across each feature. (2) For each example, Euclidian distance between the newly assigned mean feature values and the original feature values is computed. (3) Summing the individual deviations and averaging over all examples gives the total reconstruction error. The computed reconstruction error changes dynamically as new examples are processed. • Dependence on a sequence makes the algorithm as weak learner. To resolve this, probability distribution of reconstruction errors generated from a sample of random sequences is estimated. Empirically, it is found that with a maximum probability value, there exists a low reconstruction error value to which the algorithm converges. This low reconstruction error points to the adequate number of clusters also. Lastly, the algorithm currently does not work on the merging of the clusters. Section 2 discusses the algorithm in detail followed by implications and analysis of convergence in section 3. Section 4 discusses empirical results followed by conclusion.

Transductive-Inductive Learning Algorithm

Let the examples (z i = x i ) be sampled from the example set randomly without repetition in a sequential manner, thus forming a sequence Z i . The algorithm works in alternative steps by (a) creating new clusters using 1-Nearest Neighbour (NN) transductive learning and (b) learning the association of a new example to an existing cluster by evaluating the MCI (i.e inductive learning). In case no associations are found to exist, a new cluster is created by employing 1-Nearest Neighbour (NN) transductive learning. Thus the algorithm starts with no clusters at all and uses NN to initialize a new cluster. Once a cluster is initialized (say with x i and x j ), the size of the cluster depends on the number of examples getting associated with it. The MCI controls the degree of uniformity of different feature values of examples that constitute the cluster. The association of the example to a cluster happens as follows: Let the new random example (say x t ) be considered for checking the association to a cluster. If the spread of example x t from E q (x) (the mean of the q th cluster {x i , x j } and x is set of all examples in q th cluster), factored by the covariance matrix Σ q , is below C p , then x t is considered as a part of the cluster. Using MCI, it boils down to:

P{(x t -E q [x i , x j ]) T Σ -1 q (x t - E q [x i , x j ]) ≥ C p } ≤ N Cp or P{(x t -E q [x i , x j ]) T Σ -1 q (x t -E q [x i , x j ]) < C p } ≥ 1-N Cp
Satisfaction of this criterion suggests a possible cluster to which x t could be associated. This test is conducted for all the existing clusters. If there are more than one cluster to which x t can be associated, then the cluster which shows the minimum deviation from the new random point is chosen. Once the cluster is chosen, its size is extended to one more example i.e. x t . The cluster now constitutes {x i , x j , x t } and its E q (x) and Σ q recomputed. In case of failure to find any association, the algorithm employs 1-NN transductive algorithm to find a closest neighbour (in unseen data), of the current example under consideration. This neighbour together with the current example forms a new cluster. Thus the step of cluster formation or association is repeated online to reconstruct the information content in the dataset.

Implications

In this paper, the term decomposition is synonymous to cluster. The proposed work uses MCI and the probability associated with it to define a decomposition as follows:

Definition 1. Let the q th cluster D q be a decomposition, then:

D q = {x i |∀i ((x i -E q (x)) T Σ -1 q (x i -E q (x))) < C p }
A decomposition expands by testing a new point x t via the inequality (x t -E q (x)) T Σ -1 q (x t -E q (x)) < C p . The probability of the satisfaction of this inequality is given by P{(x t -E q (x)

) T Σ -1 q (x t -E q (x)) < C p } ≥ 1 -N
Cp . Thus at any point in time, while processing the dataset online, the probability of existence of D q is lower bounded by a value of 1 -(N /C p ). Lemma 1. For any cluster q, the probability of existence of D q is lower bounded by a value of 1 -(N /C p ).

It is important to note that C p plays a major role in deciding whether a new example x t belongs to D q . If the new example cannot be associated to a particular decomposition, then it is tested with other decompositions using MCI with same value of C p . Thus C p acts as a constraint in checking the homogeneity of D q . Now, since all examples (except for those used to from new clusters where the NN is used) are tested using MCI and finally associated with one or the other decomposition for which the constraint is satisfied, the total number of generated clusters is limited. This limitation is enforced indirectly via the C p . In this scenario, it can be expected that the total number of clusters is upper bounded by M/C p (M being the total number of examples in the dataset). Since the existence of a homogeneous decomposition is bounded probabilistically, the reconstruction error associated with all the examples present in the decomposition are also bounded. For all decompositions, the summation of reconstruction errors is also bounded. Thus for a particular value of C p a proof of convergence is needed for the error rates as the number of processed examples and the number of clusters increase.

Two error rates are computed as the random sequence of examples get processed. Let x i ∈ R N (N is the number of features) be in the dataset. Since the example is assigned to a particular decomposition D q , it gets a value of the mean of the all examples that constitute the decomposition. Thus the reconstruction error for the example turns out to be ||x i -E q (x)|| 2 . For each cluster q, the reconstruction error is Err Dq = n i=1 ||x i -E q (x)|| 2 (n is the number of examples in the q th cluster). As new examples are processed based on the conceptualized information from the previous examples, the total error computed at after processing the first pt cntr examples in a random sequence is

Err val = clustercntr q=1
Err Dq (cluster cntr is the total number of clusters generated after the pt cntr examples have been processed). The error rate for these pt cntr examples is Err 1 = Err val /pt cntr . Finally, the rate of error after every new cluster formation is also computed. This error is denoted by Err 2 (i.e. Err 2 = Err val /cluster cntr ).

Theorem 1. From a dataset, if examples are selected randomly and processed online in a sequential manner without repetition for a particular value of C p using the TI, then the reconstruction error rate Err 1 converges asymptotically with a probabilistically lower bound or confidence level of 1 -N /C p or greater.

Proof. Since C p defines level of compromise in information content via lemma 2 and the decompositions D q is almost homogeneous, all examples that constitute a decomposition have similar feature values. Due to this similarity between the feature values, the non-diagonal elements of the covariance matrix in the inequality above approach to zero or smaller values. Thus, Σ -1 q approaches a diagonal matrix. Multiplying a vector with diagonal matrix scales the vector by some constant factor (say S q ). Thus, if Σ -1 q ≈ S q × I, were I is the identity matrix then the inequality equates to: (x t -E q (x)) T S q I(x t -E q (x))

C p ⇒ (x t -E q (x)) T I(x t -E q (x)) Cp Sq ⇒ ||x t -E q (x))|| 2 Cp Sq
Thus, if x i = x t was the last example to be associated to a decomposition, the reconstruction error ||x i -E q (x)|| 2 for that example would be upper bounded be q , but in the worst case scenario, if the decomposition with the smallest constant factor is substituted for every other decompositions, then the upper bound on the error is Cp×Σ cluster cntr q=1

Σ n i=1 (1) S smallest ×ptcntr which equates to Cp S smallest . Since ∀q, D q s are probabilistically bounded by the Chebyshev inequality, the error rate converges with the probability associated with the Chebyshev inequality. Q.E.D

Empirical Results

The converging error rate Err 1 is depicted in left of figure 1 for a 64×64 patch in [START_REF] Sinha | Bounded multivariate surfaces on monovariate internal functions[END_REF] of starfish image. The error rate Err 2 is the computation of error after each new cluster is formed. The proof for upper bound on Err 2 follows on similar lines. For C p = 7, the error rate Err 2 is depicted on right in figure 1 for same common patch. Intuitively, it can be seen that both the reconstruction error rates converge to an approximately similar value.

Hitherto, reconstruction error and the number of clusters is dependent on a sequence presented to the learner. This points to the problem of whether an image can be reconstructed at a particular C p where there is a high probability of finding a low reconstruction error and the number of clusters, from a sample of sequences. The existence of such a probability value would require the knowledge of the probability distribution of the reconstruction error over increasing [START_REF] Vovk | Algorithmic learning in a random world[END_REF] number of examples and (2) number of clusters generated. KDE is used to estimate the probability distribution of the reconstruction error Err 1 and Err 2 . The KDE is based on a normal kernel using a window parameter that is a function of the number of points. The density is evaluated at 100 equally spaced points that cover the range of the data. The KDE empirically point to the least error rates with high probability. It was found that the error rates Err 1 , Err 2 and the number of clusters, all converge to a particular value (33.1762, 35.9339 and 38, respectively), for a given image (figure 2).

Conclusion

Given a random sequence of examples which are processed in an online sequential manner, it is possible to converge on a reconstruction of the information content of the whole dataset to a certain degree compromise with low reconstruction error using the proposed algorithm.

Lemma 2 .

 2 The value of C p reduces the initial representation of information content present in M examples to a representation of information content present in the decompositions whose number is upper bounded by M/C p .

  CpSq . Consequently, the total error after processing pt cntr examples is also upper bounded, i.e. Err val = the error rate Err 1 = Err val /pt cntr is also upper bounded. Different decompositions may have different Σ -1

Fig. 1 .

 1 Fig.1. Error rate on left(right) -Err1(Err2) for a particular sequence representing a 64 × 64 patch of starfish image[START_REF] Sinha | Bounded multivariate surfaces on monovariate internal functions[END_REF] with increasing number of processed examples(clusters) with Cp = 7.

Fig. 2 .

 2 Fig. 2. The probability density estimates for (a) Err1 (b) Err2 and (c) the number of clusters generated over 1000 random sequences representing a 64 × 64 patch of starfish image [4] with Cp = 10.

Acknowledgement

The author thanks Dmitry Adamskiy Phd candidate and Dr. Ilia Nouretdinov research assistant at the CLRC, Royal Holloway University of London and Dr. David M. J. Tax at the PRB Lab, Delft University of Technology, The Netherlands and the anonymous reviewers for their valuable and critical feedbacks on improving the manuscript.


Tatiana Tambouratzis
email: tatianatambouratzis@gmail.com

Marianna Tzormpatzaki

An Implementation of the Digital Music Stand for Custom-Made On-Screen Music Manuscript Viewing

Keywords: digital music stand, optical music recognition, music manuscript, on-screen presentation, custom-made viewing, minimal processing

An implementation of the digital music stand (DMS) is presented in this invited tutorial. The DMS constitutes an optical music recognition (OMR) system whose principal aim is the custom-made on-screen presentation of music manuscripts (MMs). Aiming at performing the minimum necessary amount of MM processing, the DMS constitutes middle ground between full-processing classical commercial OMR systems/platforms and pure MM presentation systems. The DMS is intended for musicians and choir singers at all levels of music erudition, in both cases including those with visual impairments which cannot be corrected with glasses or contact lenses; music activities ranging from study, rehearsal and performance of individuals as well as of orchestras are accommodated. The DMS implementation stages are detailed, with performance accuracy, which is necessary for correct MM presentation, also being reported.

Introduction

In the more than 50 years of optical music recognition (OMR) research, a variety of applications have been developed including (a) the preservation of old music manuscripts (MM) editions; (b) their re-editing and translation into a desired font, type of notation (for the same or a different musical instrument) and key; (c) their conversion into WAV format; (d) their classification and compilation into music databases. All of these applications aim at helping music students and researchers, performers and conductors during music study, rehearsal and performance. The enormous volume of research performed (e.g. [START_REF] Pruslin | Automatic Recognition of Sheet Music[END_REF][START_REF] Prerau | Computer Pattern Recognition of Standard Engraved Music Notation[END_REF][START_REF] Prerau | Do-Re-Mi: a Program that Recognizes Music Notation[END_REF][START_REF] Blostein | A Critical Survey of Music Image Analysis[END_REF][START_REF] Reed | Automatic Computer Recognition of Printed Music[END_REF][START_REF] Bellini | Optical Music Sheet Segmentation[END_REF][START_REF] Droettboom | Optical Music Interpretation[END_REF][START_REF] Rossant | A Global Method for Music Symbol Recognition in Typeset Music Sheets[END_REF][START_REF] Bainbridge | A Music Notation Construction Engine for Optical Music Recognition[END_REF][START_REF] Rossant | Robust and Adaptive OMR System, Including Fuzzy Modeling, Fusion of Musical Rules, and Possible Error Detection[END_REF][START_REF] Szwoch | Guido: a Musical Score Recognition System[END_REF], for a review the reader is referred to [START_REF] Rebelo | Optical Recognition of Music Symbols: a Comparative Study[END_REF]) has resulted in (a) commercially available systems such as MidiScan (http://www.musitek.com/); SharpEye and PhotoScore (http://www.visiv. co.uk/, http://www.sibelius.com/products/photoscore/ultimate.html and http://www.sibelius.com/products/photoscore/ lite.html); forScore (http://www.forscoreacom/) and platforms, e.g. Gamera (http://gamera. sourceforge.net/), used for creating born-digital MMs via the direct input of music symbols by the user (b) notation languages for recognising each music symbol of a -usually borndigital -MM and for, subsequently, representing the symbol by the corresponding prototype symbol as well as its relationship to other music symbols of the MM, for instance ABC [9-10] (c) pure MM presentation applications for on-screen viewing of MM pages, e.g. forScore (http://www.forscoreapp.com/). The digital music stand (DMS) was introduced in 2008 [13-14] as an OMR system focusing upon the custom-made on-screen and on-line reconfigurable MM presentation. This aim is accomplished via the minimal possible amount of image (MM page) processing, thus rendering the DMS a representative of the middle ground that lies between full-processing classical commercial systems/platforms and pure MM presentation applications. The DMS is intended for musicians at all levels of music erudition as well as at choir singers, in both cases including those with visual impairments which cannot be corrected with glasses or contact lenses; music activities ranging from study, rehearsal and performance of individuals as well as of orchestras are accommodated.

A complete DMS implementation is presented in this invited tutorial, which is organised as follows: Section 2 introduces the main characteristics of the DMS, contrasting them with those of existing OMR implementations; the key stages of the DMS are described in detail in Section 3, with the performance accuracy of each stage also being reported; finally, Section 4 concludes the paper.

The Main Characteristics of the DMS -Comparisons with Existing OMR Implementations

Unlike classical OMR implementations ((a)-(b) of Section 1), but like pure MM presentation applications ((c) of Section 1), no prior MM alignment, noise and/or staff line removal, or full music symbol recognition are performed. Instead, only the clefs, global accidentals and time signatures are recognized so as to allow correct on-screen MM viewing, independent of the desired MM on-screen layout and current viewing conditions.

The DMS processes scanned as well as born-digital MMs in the same uniform manner, independent of the image dimensions, format and quality. Thus, MM pages that contain noise (white or coloured) and other imperfections that are due to the state of the original MM and/or the scanning process (e.g. badly aligned, or derived from damaged MMs) are handled exactly as are those derived from crisp, near-perfect born-digital MMs.

Three viewing conditions are available to the user for custom-made on-screen MM presentation:

• Static viewing. The user specifies the desired size of the music notation, whereby the MM is automatically arranged into screens of appropriately-configured sequences of measures of the MM. The user can navigate to next/previous screens by pressing pre-specified keys; additionally, on-line changes to the music notation size are possible. • Dynamic viewing. The user specifies the desired size of the music notation as well as the tempo at which the MM is to be played. Further to static viewing, the internal metronome of the DMS causes automatic page turning as soon as the total time required for playing the part of the MM appearing on-screen has elapsed.

• Scrolling, where -as for dynamic viewing -the user specifies the desired size of the music notation and the tempo at which the MM is to be played. In this case, however, the MM is presented in a single-system format, with consecutive measures and systems sliding in a right-to-left direction and at a constant speed that is determined by the tempo. By minimizing head and eye movement, scrolling is a great aid to visually impaired musicians who can keep their fixation point constant on the centre of the screen while viewing the MM.

The DMS Key Stages

The correct operation of the DMS is centred upon the maximally accurate detection of the staff lines (stage 1). Next, each MM page is segmented independently into its constituent staves and grand staves (stage 2), as well as into the measures that appear within each of them (stage 3). The key music symbols (clefs, global accidentals and time signatures) of the MM are subsequently isolated and identified (stage 4). Finally, the key music symbols and the remaining parts of the MM are put together in such a manner as to allow custom-made on-screen viewing of the MM that can be altered online upon request (stage 5).

The reported results that are shown at the end of each sub-section have been derived from a total of 137 MMs, 115 of which are scanned and the remaining 22 are born-digital.

Staff Line Detection and (Grand) Stave Segmentation -Stages 1-2

These stages constitute the most crucial stages of the DMS implementation and are, thus, presented together. The following steps are performed:

• An initial estimation of the staff lines as well as of the staff line distance (SLD) is provided by the horizontal projection of the MM page. • Missing staff lines are filled in for each staff which has less than five detected staff lines; the SLD is used to this end. • Any extraneous staff lines are deleted from each staff which has more than five detected staff lines; again, the SLD is employed. • Non-detected staves are identified, based on the free available area as well as on the proportion and spacing of black horizontal line segments appearing within the area.

• The staff line positions are optimised within each identified staff, using the SLD as well as the location of black horizontal line segments of the original image. • Any staves that were detected during the previous step but do not correspond to real staves (e.g. due to the existence of significantly broken lines, or of unacceptable fluctuations in distance between consecutive staff lines of the same staff) are removed. Staff line and staff detection accuracy amounts to 95.95% and 90.5% for the scanned MMs and to 97.8% and 100% for the born-digital MMs.

Measure Segmentation -Stage 3

The measures within each (grand) staff are identified using the vertical projections of windows containing the grand staff and each of the two constituent staves. Measure detection exceeds 94% for both kinds of MMs.

Key Music Symbol Identification -Stage 4

This stage is vital for the accurate on-screen presentation of the MM in that it makes sure that the current key music symbols (clef, global accidentals and time signature1 ) pertaining to the viewed part of the MM are presented to the user, thus allowing him/her correctly study and/or play the MM.

A database of each key music symbols has been compiled from the available MMs. Between eight and 13 exemplars from each key music symbol have been employed to this end, each exemplar covering different sizes and fonts. Once a MM is input in the DMS, windows (of size equalling each exemplar) of the measures that were isolated in the previous stage are exhaustively matched against each exemplar and the match locations as well as the identity of the matched exemplar are retained.

Clef, global accidental and time signature identification accuracy reaches 99.8%, 92.3% and 100%, respectively, for both kinds of MMs.

MM Presentation -Stage 5

The sequence of measures of the MM is presented on-screen according to the user requirements concerning notation size, tempo and viewing condition. The key music symbols are positioned at the left end of each staff, with their identity changing whenever there is a change of clef, global accidental and/or time signature in the original MM.

DMS Evaluation

The DMS is centred upon the maximally accurate detection of the staff lines (stage 1). Subsequently, each MM page is independently segmented into its constituent staves and grand staves (stage 2) as well as into the measures that appear within each (grnad) staff (stage 3). Next, the key music symbols (clefs, global accidentals and time signatures) of the MM are isolated and identified (stage 4), thus allowing the concatenation of the appropriate key music symbols and the remaining parts of the MM in such a manner as to allow custom-made on-screen viewing of the MM (stage 5). A number of demontrations of the constituent stages as well as of the operation of the DMS will be presented at the Eight AIAI Conference.

The key music symbols affect the way the MM is played in terms of tonality as well as tempo.


Evangelia I Zacharaki

Guray Erus

Anastasios Bezerianos

Christos Davatzikos

Fuzzy multi-channel clustering with individualized spatial priors for segmenting brain lesions and infarcts

Keywords: brain tissue segmentation, lesions, infarcts, fuzzy clustering, outlier detection, MRI

Quantitative analysis of brain lesions and ischemic infarcts is becoming very important due to their association with cardiovascular disease and normal aging. In this paper, we present a semi-supervised segmentation methodology that detects and classifies cerebrovascular disease in multi-channel magnetic resonance (MR) images. The method combines intensity based fuzzy c-means (FCM) segmentation with spatial probability maps calculated from a normative set of images from healthy individuals. Unlike common FCM-based methods which segment only healthy tissue, we have extended the fuzzy segmentation to include patient-specific spatial priors for both pathological conditions (lesions and infarcts). These priors are calculated by estimating the statistical voxel-wise variation of the healthy anatomy, and identifying abnormalities as deviations from normality. False detection is reduced by knowledge-based rules. Assessment on a population of 47 patients from different imaging sites illustrates the potential of the proposed method in segmenting both hyperintense lesions and necrotic infarcts.

Introduction

Brain lesions might be associated with numerous diseases, such as multiple sclerosis (MS), traumatic brain injury, or several cardiovascular symptoms, such as arterial fibrillation, impaired cognition and others [START_REF] De Groot | Cerebral white matter lesions and cognitive function: the Rotterdam Scan Study[END_REF] [START_REF] Rovaris | Relevance of Hypointense Lesions on Fast Fluid-Attenuated Inversion Recovery MR Images as a Marker of Disease Severity in Cases of Multiple Sclerosis[END_REF]. Many methods in the literature have been developed for the automatic segmentation of brain lesions, especially MS lesions, using MR imaging [3][4]. Supervised segmentation strategies use features extracted from manual delineations as a training set to build the segmentation model. Those strategies require substantial human effort and are prone to intra-and interobserver variability. Semi-supervised strategies model the normal data and then detect lesions as outliers that are not well explained by the model [START_REF] Zacharaki | Abnormality segmentation in brain images via distributed estimation[END_REF]. These methods usually derive a prior classification from a digital brain atlas that contains information about the expected location of white matter (WM), gray matter (GM) and cerebrospinal fluid (CSF) [START_REF] Van Leemput | Automated segmentation of multiple sclerosis lesions by model outlier detection[END_REF]. The available digital atlases, such the SPM [START_REF]Statistical Parametric Mapping[END_REF] or ICBM [START_REF]International Consortium for Brain Mapping[END_REF] tissue probability maps, have very low resolution, thus are suitable only for the detection of large lesions.

Although there have been quite many attempts to segment white matter or MS lesions, relatively less attention has been given to the segmentation of stroke related infarct lesions in elderly individuals, and Alzheimer's disease or diabetic patients. Cortical infarcts are characterized by a necrotic part with intensity similar to CSF imposing ambiguity in segmentation. Shen et al. [START_REF] Shen | An improved lesion detection approach based on similarity measurement between fuzzy intensity segmentation and spatial probability maps[END_REF] applied a modified Fuzzy cmeans (FCM) clustering algorithm to classify voxels in an image according to intensity and prior probability of healthy tissue. They subsequently segment abnormal areas by using similarity to measure residual difference between the calculated tissue probabilities maps and the original prior probabilities. In such an approach lesions will be originally labeled as CSF or GM and will mislead the computation of class parameters.

Similarly to [3][9] we propose to segment infarcts within an outlier detection framework where the inconsistency between image intensity and spatial location is modeled. The proposed framework however differs in that it incorporates individualized spatial priors for both healthy and pathological tissue. The spatial priors for the healthy tissues are calculated by exploiting statistics from a healthy population group, whereas the spatial priors for the abnormal tissues are constructed in an outlier detection scheme [START_REF] Lladó | Segmentation of multiple sclerosis lesions in brain MRI: A review of automated approaches[END_REF]. For this purpose, all images from the healthy population group and testing subjects are firstly registered in a common space. The method furthermore classifies the abnormal tissue into gliotic (hyperintense) and necrotic (hypointense) regions based on multi-channel MRI (T1w and FLAIR). The five tissue types are obtained from FCM clusters by maximum membership segmentation followed by the application of some knowledge-based rules which reduce false detections.

Although the focus of the paper is to detect stroke-related infarcts and separate the hyperintense part of the infarct lesion (IL) from the hypointense infarct necrosis (IN), the method is also applied for the segmentation of other pathologies with high FLAIR signal appearance, such as white matter lesions, and these are marked also as IL.

The method has been assessed on 47 brain MR images with lesions and compared against expert manual segmentations.

Methods

Brain abnormalities manifest themselves due to either unusual intensity characteristics (evident in appropriate image modalities) or due to inconsistency between tissue intensity and spatial location. Under this scope the method applies a modified FCM algorithm which combines intensity based similarity with prior anatomical knowledge. Specifically, image segmentation is performed in two steps. First the healthy and pathological regions are extracted by applying a semi-supervised approach using healthy images as training data. Then the segmented regions are smoothed and used as patient-specific prior tissue probability maps in a modified FCM algorithm. The first step is described in section 2.1, whereas the modified FCM algorithm is presented in section 2.2.

2.1

Construction of patient-specific spatial priors

The methodology for calculating prior tissue probability at each voxel is illustrated in Fig. 1 and can be separated into two branches, the construction of spatial priors for normal classes and abnormal classes.

The spatial priors for abnormal classes are constructed following an outlier detection approach. Assuming that voxel intensities are statistically independent and follow a multivariate normal distribution we build a statistical atlas from a training set of healthy subjects. Let be the number of parameters used (MR sequences). For each voxel the statistical model consists of a 1 vector with the mean gray level intensities and a covariance matrix. Abnormalities are then detected as deviations from this model using as distance metric the voxel-wise Mahalanobis distance. This distance map, expressing probability of each voxel to be outlier, is then thresholded by a spatially varying threshold. The threshold is determined for each voxel by the ability of the statistical model to explain the training data. Specifically the average voxel-wise Mahalanobis distance is calculated for the healthy population group. Let's note it

. The applied voxel-wise threshold is • , where is a scalar (e.g.

3). This methodology extracts all abnormal regions including ILs, and INs. The abnormal regions are subsequently classified into ILs and INs based on intensity characteristics. Voxels with FLAIR intensities higher than the average intensity are classified as ILs, whereas hypointense voxels are classified as INs.

K-means segmentation WM, GM, CSF IL, IN

Smooth and merge

Spatial priors

The healthy tissue priors are calculated by segmenting the patient's T1w image into WM, GM and CSF using the k-means algorithm. Since abnormalities are present, the algorithm is not randomly initialized because abnormal tissue might influence the cluster center estimation of healthy tissue. The estimation is thus guided by the healthy training data. Specifically, the healthy images are segmented by k-means, one by one, updating each time the estimate of the clusters centers location. The average (over all subjects) center location for WM, GM and CSF clusters is finally used as initial estimate in the segmentation process of the patient's image. Therefore, the effect of the presence of abnormalities is reduced.

All crisp tissue segmentations are converted to prior probabilities by gaussian smoothing. If Ω is the set of voxel locations in the image domain and the number of classes, e.g.

5, all probability maps are collected together and normalized to satisfy ∑ 1 , for each voxel ∈ Ω. If a class is missing and not detected in the outlier detection step described above, e.g. necrosis might not be present in every image, is correspondingly reduced. These maps are utilized as individualized spatial priors in equation ( 1) as shown in the next section.

Modified FCM for healthy and abnormal tissue segmentation

The FCM algorithm [START_REF] Bezdek | Pattern recognition with fuzzy object function algorithms[END_REF] is based on minimizing an objective function with respect to a fuzzy membership U and set of cluster centroids V. Let ∈ , be the feature representation (e.g. multi-channel intensities with

2) of voxel ∈ Ω. FCM clusters the data by computing the fuzzy membership  at each voxel to the i-th class, such that ∈ 0,1 and ∑ 1, ∀ ∈ Ω. The number of classes C is assumed to be known. In the conventional approach the objective function expresses the degree of intensity based dissimilarity between the data values and the cluster centroids,  . In this modified approach the objective function includes an additional term reflecting the distance between the fuzzy membership and prior tissue probability, p ij at each voxel , which is calculated as described in the previous section. Thus is expressed as follows:

, ∑ ∑ , ∈Ω (1)
where is a distance function, m ∈ (1,∞) is a weighting exponent on each fuzzy membership controlling the degree of fuzziness and , 1 … , is a normalization weight controlling the contribution of each distance term per cluster. The use of different weights allows controlling the penalty on dissimilarity for each class. For unbiased results can take the default value of 1, as performed in this study. For the two terms in the objective function to be comparable (intensity dissimilarity and spatial inconsistency) the intensity values should be scaled in the range 0,1 , similarly to the probability range.

The constrained optimization of is expressed using Lagrange multipliers and and are determined by setting the derivative of to zero. For =2 (as usually chosen in image segmentation), the fuzzy membership and the cluster center are calculated as follows:

∑ ∑ , (2)
∑ ∑ , (3) where ,
. The algorithm iteratively optimizes the objective function in

Eq. ( 1) with the continuous update of and , until convergence is reached. The prior tissue probability maps are utilized as the membership at first iteration. In regular FCM, is the Euclidean distance. In order not to restrict the detection to data classes with hyper-spherical shape, we used the Mahalanobis distance expressed as follows:

, ′ • • (4) where is the inverse covariance matrix of the total dataset in this application. An improved implementation (to be investigated in the future) would include a different covariance matrix for each class and thus require iterative updates of based on the class membership at each iteration [START_REF] Liu | Fuzzy C-Means Algorithm Based on Standard Mahalanobis Distances[END_REF].

Knowledge-based rules

Knowledge-based rules are finally applied as post-processing steps to automatically reduce false positives. Abnormal voxels that are well explained by the statistical model of healthy anatomy, i.e. exhibiting small Mahalanobis distance, are assumed to be false detection and are eliminated. Single ischemic voxels are also removed. Two rules are applied for necrosis. Small necrotic clusters that lie across the GM-CSF interface are expected to be GM or CSF wrongly detected as IN due to partial volume effect. Necrosis is also constraint to be in proximity with a gliotic hyperintense region (IL). In contrast to other lesion segmentation methods [START_REF] Van Leemput | Automated segmentation of multiple sclerosis lesions by model outlier detection[END_REF] we didn't restrict the ILs to be situated inside WM in order not to further reduce generality.

Data description and preprocessing

The data consist of transaxial T1w and FLAIR scans obtained from different imaging sites. The acquisition of FLAIR (template space) was performed with a 3-mm slice thickness, no slice gap, voxel size=0.94×0.94×3 mm 3 , and in plane matrix 256×256. The images which didn't show any visible pathologies (72 in total) were used as the training data for statistical modeling of the healthy tissue. Different sets of data were selected for testing the method including lesions with or without necrosis. The first set (dataset1) consists of 35 diabetic patients with ILs (in 33 subjects) and necrotic infarcts (in 2 subjects). The second set (dataset2) involves 12 subjects selected from an ageing population (all subjects have both ILs and INs). The abnormality masks used for assessing the method were manually delineated by a single rater for dataset1 and by two independent raters for dataset2. Hyperintensities in FLAIR images with respect to surrounding healthy white matter tissue were marked as IL. Cortical infarcts were separated into two regions. Their necrotic part was segmented as IN, while the hyperintense rim, which usually surrounds the necrotic part, was marked as IL.

All images (from the testing subject and the normal subjects used for building the statistical atlas) were preprocessed in the same way by applying skull stripping, bias correction, global histogram matching and registration. For each subject the T1w image was first registered to the FLAIR image by affine registration. Then, for comparison of brains across different subjects, all subjects were spatially normalized to one of the subjects chosen as a template, by registering their T1w images with the template's T1w image using the HAMMER registration algorithm [START_REF] Shen | HAMMER: hierarchical attribute matching mechanism for elastic registration[END_REF]. HAMMER algorithm achieves registration by hierarchically finding correspondences on edge voxels with high (feature-based) similarity. Voxel pairs with low similarity, such as the ones around pathology, are not included in the deformation process. Therefore, although the algorithm is designed for images without pathology, it has been observed that it is also quite robust when the images include small localized abnormalities. The FLAIR images and the expert-defined masks were subsequently warped to the template space using the deformation fields estimated for their corresponding T1 images.

Assessment

The results of the automated method were compared with lesion and infarct delineations by human experts. The segmentation of IL was assessed using the Dice coefficient as measure of spatial overlap. The Dice coefficient over two segmentations, S 1 and S 2 , is the volume of intersection of the two segmentations divided by the mean of the two segmentation volumes:

2 • | ∩ |/ | | | | .
The necrotic part of infarcts is usually too small to be appropriately assessed by criteria measuring volume overlap, such as the Dice coefficient. Thus necrotic voxels were grouped into clusters based on connectivity and segmentation assessment was performed based on the number of correctly classified clusters. A necrotic cluster is considered as a correct detection if the voxels within the cluster partially overlap with manually segmented necrosis. Clusters with a volume smaller than a threshold Vmin were not considered in the evaluation process.

Results

Segmentation of lesions and infarcts

Fig. 2 (last row) illustrates some examples of IL and IN segmentation by the proposed method. The results are compared against manual delineations (IL available for subjects 3 and 4, and IN available for subjects 1, 2, 3). In case of large lesions overall good agreement is observed. In order to quantitatively assess lesion segmentation, Dice coefficients have been calculated for both datasets and plotted against the total lesion load (or average volume in the case of two raters) for all subjects in Fig. 3. The inter-rater variability is also calculated for dataset2 by defining the one of the two masks as ground truth and assessing the other. It can be observed that the agreement between automatic and manual segmentation increases with increasing lesion load. This is expected since small lesion clusters might be missed by the method, due to partial volume effects, smoothing and/or warping of the images and manual masks. We should note that any further volumetric measurements for clinical assessment should be performed in the original space by warping back the images with segmented lesions. The segmentation of infarct necrosis is assessed based on the number of correctly or incorrectly detected clusters, defined as true positives (TP) and false positives (FP), over all subjects. The effect of the cluster volume threshold, V min , on segmentation accuracy, is shown in Fig. 4. On dataset1 (Fig. 4a) seven necrotic clusters were in total manually outlined for the two subjects. Two of them are missed by the automated method (for V min =0). The volume of each of the missed clusters is equal to a single voxel; thus it is reasonable to miss these clusters. The assessment of detection rate is different for dataset2 where two manual masks (R1 and R2) are available. First we assumed R1 is ground truth and calculated the sensitivity (TP rate) of the automated method as the ratio of correctly detected clusters by the method versus the number of clusters correctly detected by the 2 nd rater (after volume thresholding). Then we assumed R2 is ground truth and compared our method against the accuracy of the 1 st rater. The sensitivity is plotted in Fig. 4b. The results show that the method performs more similar to rater 2 in necrosis segmentation, and more similar to rater 1 in IL segmentation (Fig. 3

Discussion and Conclusions

We presented a method that uses multi-channel MR images to segment healthy and abnormal brain tissue employing a statistical model, together with a modified FCM algorithm. The segmentation of enhancing lesions approached the accuracy of manual segmentation for high lesion load. The detection of infarcts was only slightly worse than the performance of human experts. However more validation is required because there was a considerable disagreement in spatial correspondence between expert segmentations. The assessment of the method on synthetic databases might provide useful insight. Also the sensitivity of the method to the selected parameters should be investigated.

The method was applied to segment and classify vascular pathology; it can nevertheless also be applied to detect lesions that have not vascular origin. Adaptation of the knowledge-based rules might be required in this case.

Extensions of the FCM algorithm have been proposed in the literature to reduce the algorithm's sensitivity to noise and intensity inhomogeneity. The extensions incorporate spatial regularisation, topological constraints or bias field correction, as briefly reviewed in [START_REF] Ji | A modified possibilistic fuzzy c-means clustering algorithm for bias field estimation and segmentation of brain MR image[END_REF]. In this paper we did not apply any of those extensions since we focused on the lesion and infarct segmentation. It is expected however that the improved modeling of the healthy tissue will improve also the segmentation of the abnormal tissue.

Fig. 1 .

 1 Fig. 1. Methodology of the calculation of spatial priors for normal tissue (white matter (WM), gray matter (GM) and cerebrospinal fluid (CSF)), and abnormal tissue (hyperintense infarct lesion (IL) and infarct necrosis (IN))

Fig. 2 .

 2 Fig. 2. Visual segmentation assessment on 4 subjects. From top to bottom: FLAIR with manually outlined IL (available only for subjects 3 and 4), T1w with manually outlined IN, automated segmentation (blue: CSF, dark green: GM, light green: WM, pink: IN, white: IL). The manual segmentations are shown in color scale with darker color indicating the presence of IL or IN in neighboring slices.

Fig. 3 .Fig. 4 .

 34 Fig.3. Dice scores in relation to total IL volume (in voxels) for each subject in dataset1 (left) and dataset2 (right). Comparison of our method (S) with manual segmentations by raters R1 and R2, and inter-rater variability (R1 vs R2) are shown for dataset2. In this case the total volume represents the average volume.

  , right). The total number of clusters (summed over all subjects) and their median volume (in voxels) are shown in Fig.4c. The statistics show that the necrotic part of the infarcts in dataset2 included many unconnected components making automated detection more difficult.

Acknowledgments

This research was supported by a Marie Curie International Reintegration Grant within the 7 th European Community Framework Programme.


Nantia Iakovidou
email: niakovid@csd.auth.gr

Eleftherios Tiakas
email: tiakas@csd.auth.gr

Konstantinos Tsichlas
email: tsichlas@csd.auth.gr

DISCO: a new algorithm for detecting 3D protein structure similarity

Keywords:

Protein structure similarity is one of the most important aims pursued by bioinformatics and structural biology, nowadays. Although quite a few similarity methods have been proposed lately, yet fresh algorithms that fulfill new preconditions are needed to serve this purpose. In this paper, we provide a new similarity measure for 3D protein structures that detects not only similar structures but also similar substructures to a query protein, supporting both multiple and pairwise comparison procedures and combining many comparison characteristics. In order to handle similarity queries we utilize efficient and effective indexing techniques such as M-trees and we provide interesting results using real, previously tested protein data sets.

Introduction

The three dimensional (3D) structure of a protein is directly related to its functionality and could probably even determine it. This is why the comparison of protein structures became such an extremely important problem in structural biology and other scientific areas as well, such as medicine, agriculture and industry. Hundreds of algorithms that perform protein structure comparison have been proposed so far, with a view to compare and group proteins by structural similarity or identify distant homologues of protein families. The main reason for using 3D comparison algorithms, instead of using for example sequence comparison algorithms, is to detect functional similarities between proteins with low sequence similarity ( [START_REF] Micheletti | MISTRAL: a tool for energy-based multiple structural alignment of proteins[END_REF], [START_REF] Lichtarge | Evolutionary predictions of binding surfaces and interactions[END_REF]).

In this paper we present a novel algorithm that performs 3D protein structure comparison at the level of C-alpha atoms and aims at detecting similarity between a specific query protein and a set of one or more known protein structures. For this purpose we use a simplified framework for representing the protein structures inspired by [START_REF] Zhi | Representing and comparing protein structures as paths in three-dimensional space[END_REF]. In our framework each protein molecule is represented as a sequence of weights, which come up from the combination of two measures. If we suppose for example that n is the number of C-alpha atoms in the protein molecule, then what we compute is a) the n-1 distances between consecutive C-alpha atoms and b) the n-2 cosines of their associated turning angles, which according to [START_REF] Zhi | Representing and comparing protein structures as paths in three-dimensional space[END_REF] consist a sufficient angle descriptor and no use of torsion angles is needed. Then we combine these two quantities n-1 times in order to produce a series of the hybrid measure w (w 1 , w 2 , ..., w n-1 ). In this way, local secondary structural information is maintained and also the use of the DIStance and COsine measures, from which our algorithm's name (DISCO) is produced, ensures that the algorithm is rotation and translation independent.

On the next step, we divide this sequence of weights into overlapping subsequences, which correspond to all possible consecutive protein substructures of the molecule. After repeating the aforementioned procedure for all proteins, we make use of indexing techniques (M-trees) in order to store these data and search for similarities either between pairs of proteins (Pair-DISCO Algorithm), either between a particular query protein and a set of other proteins (Multi-DISCO Algorithm). The use of overlapping subsequences enables us to detect not only similar structures but also similar substructures (motifs) in protein molecules, without thinking about sequentiality, which means about the order of appearance of the substructures into the molecule. This is very important because it means that we can capture both local and global similarity of protein structures and also because it enriches our algorithm with flexibility.

In fact, the use of indexing technique in our approach helps us to construct a database of protein structures in order to efficiently support and handle similarity queries of a specific protein across the totality of proteins that are stored in the determinate database. As will be shown in the paper, the indexing method calculates accurate similarities without having access to the whole database. In this way, the searching procedure becomes more efficient and rapid providing also dynamic data handling.

The rest of the paper is organised as follows. Section 2 presents previous related work in this area. In Section 3 our proposed approach and its experimental evaluation are presented. Finally some conclusions and future work are provided in Section 4.

Related Work

Numerous methods for performing protein structural searches have been proposed in recent years. Some of them perform alignment in rigid structures between amino acids at the level of C-alpha atoms. In general, all protein structure alignment programs optimize some mathematical definition of structural similarity. A typical example of this category is the root mean square deviation measure (RMSD). Although these methods provide a pragmatic definition of structural similarity, they fail to capture similar structures with extensive conformation changes such as internal rearrangements, which means that these methods ignore the flexibility of the polypeptide chains ( [START_REF] Zhi | Representing and comparing protein structures as paths in three-dimensional space[END_REF]).

Several other algorithms align proteins at the level of secondary structure elements (SSEs). Such kind of approaches require a criterion for assigning secondary structures to proteins, as belonging to helices, strands or loops or not being part of an SSE at all. Tableau-based methods generally belong to this category, such as SA Tableau Search ( [START_REF] Stivala | Fast and accurate protein substructure searching with simulated annealing and GPUs[END_REF]) and IR Tableau ( [START_REF] Zhang | A fast indexing approach for protein structure comparison[END_REF]). The truth about SSE methods is that there is not an exact procedure of assignment and also opinions vary about the precise beginning and end of SSEs ( [START_REF] Stivala | Fast and accurate protein substructure searching with simulated annealing and GPUs[END_REF]). Furthermore, using only SSEs means that regions of protein structures that are not defined as being part of an SSE are not used at all and this fact can lead to less sensitive results.

A parameter that is also taken into account during the study of 3D protein structures is sequentiality, which means subsequent amino acids in one protein must correspond to subsequent amino acids in the partner protein. The majority of methods follow this restriction while the number of methods that are nonsequential is still limited ([6], [START_REF] Kolbeck | Connectivity independent protein-structure alignment[END_REF]). The drawback of sequential approaches is that they can decrease the possibility to discover evolutionary relationships ( [START_REF] Xie | Detecting evolutionary relationships across existing fold space[END_REF]) as new protein structures can arise from the combination and permutation of substructures of a protein ( [START_REF] Fong | Modeling the evolution of protein domain architectures using maximum parsimony[END_REF]).

Methods that do not perform any kind of alignment have also been proposed, such as [START_REF] Zhi | Representing and comparing protein structures as paths in three-dimensional space[END_REF] and [START_REF] Budowski-Tal | FragBag, an accurate representation of protein structure, retrieves stuctural neighbours from the entire PDB quickly and accurately[END_REF]. These methods use various ways to represent the protein molecules for example as paths or vectors, in order to use this representation to identify candidate sets of structural neighbours. The drawback of these approaches is that they perform an approximate representation of the protein molecule and this can lead to loss of secondary structure information and less accurate results.

Certain other issues also arise when structural comparison of proteins is studied such as pairwise and multiple comparison. The result of a pairwise alignment or comparison procedure concerns two particular proteins, while the multiple comparison procedure gives similarity results between a certain protein and a list of other known proteins. Another term called 'dynamics-based-alignment' has also been recently introduced ( [START_REF] Zen | Correspondences between lowenergy modes in enzymes: dynamics-based alignment of enzymatic functional families[END_REF]) and is intended to compare the dynamic motions of different proteins. Furthermore, methods that perform structural alignment between a model protein and the true known structure of the protein have also been proposed ( [START_REF] Ortiz | MAMMOTH: an automated method for model comparison[END_REF]) and they are known as model comparison methodologies, but these last two subjects are out of the scope of the current work.

The Proposed Approach

Protein Data Representation Algorithm

Let P be a set of proteins, and let a protein P a ∈ P , and its length n = |P a |, which is defined as the total number of C-alpha atoms that it contains. For the protein P a let x[1...n], y[1...n], z[1...n] be the corresponding 3D coordinates of the C-alpha atoms that it contains in the Euclidean three dimensional space, following the order of the protein structure.

We compute the n-1 Euclidean distances between consecutive C-alpha atoms using the formula:

d i = √ (x[i] -x[i -1]) 2 + (y[i] -y[i -1]) 2 + (z[i] -z[i -1]) 2
for i = 2, ..., n. Note that d 1 = 0. At the following, we compute the n -2 convex angles θ i ∈ [0, π] between the linear vector segments u, v associated to the secondary structures of proteins, using the cosine formula:

cos θ i = u • v |u||v|
where,

u • v = (x[i -1] -x[i])(x[i + 1] -x[i]) + (y[i -1] -y[i])(y[i + 1] -y[i]) +(z[i -1] -z[i])(z[i + 1] -z[i]) |u| = √ (x[i -1] -x[i]) 2 + (y[i -1] -y[i]) 2 + (z[i -1] -z[i]) 2 |v| = √ (x[i + 1] -x[i]) 2 + (y[i + 1] -y[i]) 2 + (z[i + 1] -z[i]) 2 for i = 2, ..., n -2. Note that θ 1 = θ n-1 = 0.
We compute the convex angles between all the protein atoms sequentially according to the Protein Data Bank file format ( [START_REF] Berman | The Protein Data Bank[END_REF]). In this way we can capture the structural characteristics of all the recorded side-chains of the proteins in their regular backbone structure.

For the protein P a we normalize the calculated distances d i and angles θ i in the interval [0, 1]:

nd i = d i max i {d i } nθ i =
θ i π and we calculate the following combined linear measure for the C-alpha protein atoms:

w i = a • nd i + (1 -a) • nθ i
for i = 1, ..., n -1, and for a selected real number a ∈ [0, 1], which defines how much the distances or the angles will affect the final calculated weights w i between the protein C-alpha atoms. Finally, we construct the protein's data as a sequence of the weights w i :

P a -→ (w 1 , w 2 , ..., w n-1 )
Henceforth, all proteins of P will be represented by this sequence of weights.

Data Indexing

Among the existing indexing schemes for metric spaces, we use the Metric-Tree (M-tree) ( [START_REF] Ciaccia | M-tree: An efficient access method for similarity search in metric spaces[END_REF]) to index the protein data, as it combines both efficient query processing and dynamic data handling. The M-Tree is a balanced tree that can index objects with attributes in a metric space, compared by distance functions which satisfy the metric properties (symmetry, positivity, triangular inequality).

A simple example of such an indexing is depicted in Figure 1. Using the M-tree, we have the advantage of inserting, updating and deleting dynamically new and old protein structures, and supporting efficiently nearest neighbors and range queries. Moreover, with the defined similarity (or distance) measures, which satisfy the metric space properties, we can support similarity queries between large protein structures, by partitioning the proteins into subproteins and using incremental nearest neighbor queries for ranking. The next subsections describe the proposed methodology for such similarity queries.

Top-k Protein-to-Protein Similarity Query (Multi-DISCO Algorithm)

Let P be a set of proteins. Each protein ). We are interested in ranking all proteins in comparison with a given query protein P q , and to retrieve the best top-k results. The idea is as follows.

P i ∈ P is divided into |P i | -T + 1 sub- proteins ,
For any sub-protein S j of P q we retrieve incrementally its nearest neighbor sub-proteins. For that purpose we use an M-tree cursor, which we initialize to the first nearest neighbor of S j , and we retrieve the next nearest neighbors of S j incrementally. An M-tree cursor (which has been implemented in the M-Tree project ( [START_REF]The M-tree Project[END_REF])), provide sorted access from a specific object to the nearest neighbor order from that object and retrieve incrementally its nearest neighbors one by one. All existing sub-proteins into the M-tree will have a specific nearest neighbor order position from S j . However, we do not leave the M-tree cursor to scan the entire indexed data-set by retrieving all nearest neighbors. When we retrieve at least one sub-protein from all existing proteins, we stop scanning the nearest neighbor order with the cursor. For any retrieved sub-protein we record its nearest neighbor order position as a ranking value of its corresponding protein. We repeat this process for any existing sub-protein S j of P q and we keep the ranking results in a table. Then, the final rankings of the proteins are computed in another table R by taking the minimum or the maximum or the average of the recorded corresponding positions. The top-k similar proteins are extracted from the table R after an ascending sort. It is important to mention here that we can tune the parameter T , which is the length of each sub-protein, in order to have small or large sub-proteins.

Pairwise Similarity Query (Pair-DISCO Algorithm)

The pairwise similarity query is performed in a straightforward fashion. The idea is to take any sub-protein S1 i from the first protein P 1 and to compute its distance to any sub-protein S2 j from the second protein P 2 . Then, using a specific rule we compute the average of the minimum or maximum or average distances and return this value as a final score. When this score is close to 0, the proteins are similar. When this score is close to 1 the proteins are dissimilar. It is important to note that we retrieve the data of the proteins P 1 , P 2 , from the Mtree into memory, before the distance computations. This leads to a significant improvement in the performance of the query.

Motif Extraction (Sub-structure query)

As already mentioned before, we divide each protein P i ∈ P into |P i | -T + 1 sub-proteins, with length T . In this way, motif extraction is achieved as an intermediate step of Multi-DISCO algorithm, where we detect at least one sub-protein from all existing proteins to be the most similar to the each time examinant subprotein. Because sub-structures of proteins are overlapping, for example if the first sub-structure is [1 2 3 4 5 6 7 8 9 10] then the second one will be [2 3 4 5 6 7 8 9 10 11], we have ranking results for all possible sub-proteins of the molecule. An example is shown in Table 1. Let us consider a dataset of 20 proteins. The numbers in the first column represent the protein and its sub-protein respectively. The second column contains the most similar sub-proteins from all twenty proteins of the example dataset in a row. That is to say, the first number of the second column corresponds to the most similar sub-protein of the first molecule, the second number (which is 120) corresponds to the most similar sub-protein of the second protein and so on.

Experimental Evaluation

The best comparison results were provided using the following parameters: subprotein length equal to the 20% of the average length of the protein dataset, the weight of the hybrid measure equal to 0.05 and using the min rule to produce the final score. These numbers were extracted after extensive experimental procedures, the results of which will not be presented in the current paper due to lack of space. On the other hand, the number of proteins |P | contained in the dataset and the total number of sub-proteins N that were created after the processing of the dataset, vary according to the each time tested protein dataset.

It is already known from biology that proteins that belong to the same family (a superfamily or a sub-family) typically have similar three-dimensional structures and functions ( [START_REF] Needleman | A general method applicable to the search for similarities in the amino acid sequence of two proteins[END_REF]). We first used Pair-DISCO algorithm to compare pairs of proteins that belong to the same protein families and to the same protein subfamilies. Similarly we applied the algorithm to compare proteins that belong to different protein families. By extensive experimentation we set a threshold for Pair-DISCO algorithm to distinguish significant similar proteins from simply similar ones. In general, values near zero indicate similar proteins. Significant similarity though is detected when the value of Pair-DISCO result is under 0.05.

We then applied Pair-DISCO algorithm in a representative sequence-independent dataset ( [START_REF] Fischer | Assessing the performance of fold recognition methods by means of a comprehensive benchmark[END_REF]), obtained by using the Protein Data Bank ( [START_REF] Berman | The Protein Data Bank[END_REF]) and taking into account only structural criteria. Each structural protein family is equally represented in this dataset and also every chain within it has less than 30% sequence identity. Table 2 summarizes the results of the Pair-DISCO algorithm. The first two columns contain the PDB id of the compared proteins and the third column shows the obtained similarity score. The forth column uses the aforementioned defined threshold in order to discriminate the significant similarities.

In our results we can see that proteins with PDB ids 1onc and 7rsa, which correspond to proteins P-30 and RNase-A respectively, are not signed to be significantly similar. This happens because the first protein performs functions and activities that the second one does not possess, even though both of them are part of the same protein family (RNase) ( [START_REF] Mosimann | Refined 1.7 A X-ray crystallographic structure of P-30 protein, an amphibian ribonuclease with anti-tumor activity[END_REF]). On the other hand, proteins with PDB ids 1gal and 3cox were found to be significantly similar and this is true because they are enzymes that belong to the same category (oxidases), which constitutes a subclass of the oxidoreductases protein family. They also perform a similar function which is to catalyze an oxidation-reduction reaction, involving molecular oxygen (O 2 ) as the electron acceptor.

We also tested the performance of Multi-DISCO algorithm on numerous datasets, results of which are available, but we mainly focused on groups of proteins that were previously used as test cases ([19], [START_REF] Micheletti | MISTRAL: a tool for energy-based multiple structural alignment of proteins[END_REF]). In particular, we considered two datasets of proteins consisting of globins and serine proteinases, two of the most extensively studied protein families in biology, as shown in Table 3. Table 4 contains some representative results that derive from the application of the algorithm to the aforementioned dataset. The query protein is the one with PDB id 3est that belongs to the serine proteinases family and the ranking results using the min rule are listed on column two. Observe that the algorithm ranks first correctly the protein with PDB id 3est, which is in fact the query protein. This happens because we include every time the query proteins in the examined datasets for reasons of verification. Our algorithm also ranks correctly in the first places, proteins that belong to the same family and particularly to the same sub-family with that of the query protein. Correspondingly, the algorithm ranks to the last places the proteins that belong to the globins family.

6 Fig. 1 .

 61 Fig. 1. Data Indexing with the M-Tree.

  with length T , which are indexed in an M-tree file. Let us consider an example. If a protein P i has length 10, depicted by the series [1 2 3 4 5 6 7 8 9 10] and T = 3 then the set of |P i | -T + 1 sub-proteins that our algorithm will create will be the following: ([

Table 1 .

 1 Example of sub-structure (motif) query.

	1-1: 1 120 6 45 103 1774 182 921 62 38 50 125 29 8 12 13 7 49 20 9
	1-2: 1 80 21 4 48 1028 212 450 88 23 14 78 11 36 13 28 3 10 34 5
	1-3: 1 56 16 4 52 1016 147 423 122 34 58 43 18 27 19 20 11 10 25 3
	1-4: 1 86 3 33 43 1143 174 298 11 101 4 7 10 8 23 52 13 32 39 20
	1-5: 1 66 22 25 19 866 67 136 10 74 13 15 16 5 24 30 4 96 61 8
	1-6: 1 59 23 13 21 1166 61 258 24 117 5 16 26 4 14 17 6 15 60 20
	1-7: 1 29 24 60 46 1126 144 226 32 193 13 47 58 9 10 6 5 12 27 18
	1-8: 1 8 17 108 65 1399 96 438 66 194 48 11 146 22 12 7 5 4 21 19
	1-9: 1 4 9 51 75 1461 55 406 77 165 33 11 159 32 28 12 2 8 10 40
	1-10: 1 2 11 35 25 1081 20 303 24 157 22 7 91 46 10 5 8 9 3 21
	1-11: 1 7 4 6 23 650 8 248 18 144 63 2 38 70 20 10 11 3 13 32
	1-12: 1 7 5 3 10 473 11 347 103 194 37 4 8 34 12 16 30 2 19 55
	1-13: 1 11 6 5 9 271 15 260 104 230 47 3 10 8 14 4 36 2 24 34
	1-14: 1 11 6 5 8 263 13 282 127 339 53 4 15 10 16 2 41 3 24 44
	1-15: 1 26 2 6 29 357 58 289 63 345 16 4 10 48 24 7 19 3 21 23
	1-16: 1 78 3 59 42 468 56 376 8 282 22 6 7 37 29 13 2 12 33 27
	1-17: 1 61 3 67 37 509 26 402 6 104 22 10 7 34 14 32 2 19 15 17
	1-18: 1 84 4 75 81 330 44 446 20 95 49 16 8 40 5 9 2 12 15 23
	1-19: 1 35 2 58 141 160 57 143 41 28 27 10 38 21 8 7 3 4 77 6
	1-20: 1 22 17 69 127 858 49 160 40 14 59 15 28 11 4 3 12 2 71 10

Table 2 .

 2 Pairwise results of Pair-DISCO algorithm.

	Protein1 Protein2 Pair-DISCO Signif. Similar
	1npx	3grs	0.0198731	YES
	1onc	7rsa	0.0809089	NO
	1osa	4cpv	0.0295322	YES
	2cmd	6ldh	0.0510457	NO
	1aba	1ego	0.0567986	NO
	1eaf	4cla	0.0234494	YES
	2sga	4ptp	0.0321825	YES
	1aaj	1paz	0.0259216	YES
	5fd1	2fxb	0.0753062	NO
	1gal	3cox	0.0240907	YES
	1tlk	2rhe	0.0306053	YES
	1omf	2por	0.0508227	NO
	8i1b	4fgf	0.0645641	NO
	1mup	1rbp	0.0321935	YES
	1arb	4ptp	0.0330052	YES
	2pia	1fnr	0.0267139	YES
	3cd4	2rhe	0.0690796	NO
	2mnr	4enl	0.0190777	YES
	2gbp	2liv	0.0325059	YES
	1fxiA	1ubq	0.03846	YES

Conclusion

In this paper we proposed a methodology for detecting similarity among or between 3D protein structures using a robust indexing technique. Our proposed scheme allows both pairwise and multiple comparisons and also it is the first time that a residue-based method allows motif (substructure) queries, as well. Furthermore, the proposed algorithms are capable of capturing both global and local similarity by tuning the appropriate parameters and combine many attributes such as translations-rotations independence, non-sequentiality and flexibility.

The results of the paper constitute a first step in our research. As a future work we will include comparison of our methodology to other well-known algorithms. We also plan to investigate how our algorithms behave in a parallel system with multiple disks. Additionally, we plan to see how our algorithms could be implemented in distributed environments and in particular we aim at the Map-Reduce framework ( [START_REF] Gaggero | Parallelizing bioinformatics applications with MapReduce[END_REF]). In this way, we will be able to allow for a larger database that will maintain thousands of proteins and will support queries and updates with good response times.