Jakub M Tomczak
email: jakub.tomczak@pwr.wroc.pl

Jerzy Świątek
email: jerzy.swiatek@pwr.wroc.pl

Personalisation in Service-Oriented Systems Using Markov Chain Model and Bayesian Inference

Keywords: modeling, Markov chain, Bayesian inference, incremental learning

In the paper a personalization method using Markov model and Bayesian inference is presented. The idea is based on the hypothesis that user's choice of a new decision is influenced by the last made decision. Thus, the user's behaviour could be described by the Markov chain model. The extracted knowledge about users' behaviour is maintained in the transition matrice as probability distribution functions. An estimation of probabilities is made by applying incremental learning algorithm which allows to cope with evolving environments (e.g. preferences). At the end an empirical study is given. The proposed approach is presented on an example of students enrolling to courses. The dataset is partially based on real-life data taken from Wrocław University of Technology and includes evolving users' behaviour.

Introduction

Many modern computer networked systems (e.g. e-commerce, web services) are used to provide services to users. Moreover, if the services play a crucial role in the system, such systems could be called service-oriented systems (SOSs) [START_REF] Brzostowski | A Rough Set-Based Algorithm to Match Services Described by Functional and Non-Functional Features[END_REF], [START_REF] Grzech | Translations of Service Level Agreement in Systems Based on Service Oriented Architecture[END_REF], [START_REF] Grzech | QoS-aware infrastructure resources allocation in systems based on service-oriented architecture paradigm[END_REF]. [START_REF] Grzech | Modeling and Optimization of Complex Services in Service-Based Systems[END_REF], [8]. However, in SOSs there are several main problems that have to be concerned [START_REF] Brzostowski | A Rough Set-Based Algorithm to Match Services Described by Functional and Non-Functional Features[END_REF], [START_REF] Grzech | Modeling and Optimization of Complex Services in Service-Based Systems[END_REF]: i) user's demand formulation and contract negotiation; ii) user's demand matching with accessible services including aspects such as e.g. knowledge about users' behaviour; iii) service execution on physical machines concerning network quality of service (QoS).

To solve the mentioned problems a process consisting of negotiation, discovery, and execution stages could be proposed. The stage for demand translation and contract negotiation could be based on ontology knowledge representation [START_REF] Pastuszko | Processing and negotiation of natural language based contracts for Virtual Organizations[END_REF]. Next, in the service discovery there are three main procedures: i) service matching with user's demand (contract), ii) personalisation of services, iii) new service composition if there is no accessible service fulfilling user's demand. Service matching could be made by applying rough set theory [START_REF] Brzostowski | A Rough Set-Based Algorithm to Match Services Described by Functional and Non-Functional Features[END_REF], and service compositioni.e. using ontologies [8]. Because the personalisation is the main topic of this work, thus it will be discussed in Section 3. The last stage (execution) has to handle final aspects of the whole process and be QoS-aware [START_REF] Grzech | Translations of Service Level Agreement in Systems Based on Service Oriented Architecture[END_REF], [START_REF] Grzech | QoS-aware infrastructure resources allocation in systems based on service-oriented architecture paradigm[END_REF], [START_REF] Grzech | Modeling and Optimization of Complex Services in Service-Based Systems[END_REF]. This paper consists of following sections. In Section 2 a contribution of proposed approach to the technological innovation is presented. Next, the general problem of personalisation is described and the solution using Markov chains is outlined. Moreover, a recursive expression for decision making using Bayesian inference is presented. At the end an empirical study is conducted. The presented approach is applied to the real-life problem of students enrolment to courses at Wrocław University of Technology (WTU) and could be used as a new functionality in the existing educational platform (so called EdukacjaCL).

Contribution to sustainability

Applying the personalisation method in service-oriented systems aims in increasing quality of user service. It is widely used in the e.g. e-commerce [START_REF] Anand | Intelligent Techniques for Web Personalization, in: "Intelligent Techniques for Web Personalization[END_REF], but there is a constant need for new approaches and applications. In this paper an approach of a Markov model as a knowledge representation and Bayesian inference as a reasoning method is proposed. Furthermore, an interesting result of using Markov model and Bayesian inference is that the decision is made due to a recursive procedure (see Section 3.1). Moreover, using knowledge about users gives SOSs a new function of sustainability.

Moreover, an another novelty is using an incremental learning paradigm for probabilities estimation [START_REF] Tomczak | Bayesian Classifiers with Incremental Learning for Nonstationary Datastreams[END_REF]. Learning about users' preferences, which evolve in time, is one of the crucial aspects in modern SOSs. And to solve it some adaptive approach has to be applied. Applying incremental learning affects in sustaining model accuracy.

Therefore, in this paper a compact framework using probabilistic model and inference for personalisation problem is proposed. The framework tries to maintain the system's sustainability by making optimal decisions about users' demands and keep an up-to-date knowledge about users.

Personalisation in service-oriented systems

Mining knowledge about users in computer systems could lead to increasing quality of user service and profits to service provider(-s), e.g. to overcome so called information overload [START_REF] Anand | Intelligent Techniques for Web Personalization, in: "Intelligent Techniques for Web Personalization[END_REF], [START_REF] Pierrakos | Web Usage Mining as a Tool for Personalization: A Survey[END_REF], [START_REF] Webb | Machine Learning for User Modeling[END_REF]. Therefore, there are different aspects of personalisation, concerning e.g. recommendation, user tutoring, adaptation of layout and content, and so on [START_REF] Webb | Machine Learning for User Modeling[END_REF].

However, in this paper the main concern is put on the recommendation task. In other words, on example of an educational platform, during student's enrolment to a course, a sorted list of courses the best suited to her/his demands is presented. In the literature there are different approaches to personalisation [START_REF] Anand | Intelligent Techniques for Web Personalization, in: "Intelligent Techniques for Web Personalization[END_REF], e.g. individual or collaborative, user or item information, memory-or model-based.

In presented approach it is assumed that user's decision depends on previously made decisions. For example, a student enrols to a new course based on courses she/he has been enrolled in the past. Besides, the decision is rather uncertain because of e.g. lack of information, and that is why a new decision could be described by a probability distribution function. Hence, the Markov chain model [START_REF] Bubnicki | Pattern Recognition Algorithms for Simple Markov Chains[END_REF], [START_REF] Zhu | Using Markov Chains for Link Prediction in Adaptive Web Sites[END_REF] seems to be a proper knowledge representation to model users' behaviour because it allows to model a situation that previous decisions affects the current decision. Then, to reason about the maintained knowledge the Bayesian inference is used.

Problem statement

As it was mentioned in previous section, a users' behaviour is modelled using Markov chains. Thus, the problem of personalisation could be stated as a classification task.

Let assume that an input in the N th moment is described by a vector of features, uN ∊U=U1×U2×…×US (card{U s } = K s < 0 ℵ ). The sequence of inputs is denoted by

) ,..., , ( 2 1 N N u u u u =
. A decision in the N th time step is denoted by j N ∊M={1,2,...,M}.

Furthermore, we assume that u n , and j n , for each n = 1, 2, ..., N are realizations of stochastic processes drawn with distributions P n (j n ) i P n ( n u |j n ). However, in further considerations it is assumed that decisions j 1 , j 2 , …, j N forms a Markov chain [START_REF] Bubnicki | Pattern Recognition Algorithms for Simple Markov Chains[END_REF] described by an initial (a priori) vector of probabilities, p

( 1

i j P p i = =
, i=1,2,…,M, and a sequence of transition matrices,

P N = [p N,i,j ], i,j=1...M, ) | ( 1 , , i j j j P p N N N j i N = = = + . Besides, it is assumed that the observations of users are independent, ∏ = = N n n n n N N N j P j P 1 ) | ( ) | ( u u .
Thus, the problem could be stated as follows. Having a new user the best suited decision should be made. It could be made by minimizing following risk functional

( ) [ ] ∑ = Ψ = Ψ N n n n n N J L R N 1 ) ( , ) ( X E , (1)
where

N Ψ = [Ψ 1 Ψ 2 … Ψ N ] is a vector of decisions, ) ( n n n i u Ψ =
, E[•] is an expected value and L(•,•) is a chosen loss function. It is easy to notice [START_REF] Bubnicki | Pattern Recognition Algorithms for Simple Markov Chains[END_REF] that to minimize risk functional (1) it is enough to consider

( ) [ ] ∫ ∑ = ⋅ ⋅ Ψ = Ψ = Ψ M m n n n m n n n n n n n n m P p m L J L R 1 ) ( ) | ( )) ( , ( ) ( , ) ( u d u u X E , ( 2
)
where a priori distribution in N th moment is

∑ = - ⋅ = M l l n l m n m n p p p 1 ) ( 1 , , ) (
. Thus, we can propose an optimal decision making algorithm (Bayesian algorithm)

{ } ) , ( min arg ,... 2 , 1 n M l n l r i u = = , (3)
where:

) | ( ) , ( ) , ( 1
) ( ∑ = ⋅ ⋅ = M m n n m n n m P p l m L l r u u
, and for 0-1 loss function:

) ( ) ( ) | ( ) , ( l n n n l n n l P p l δ δ ∆ = ⋅ = u u . ( 4
)
It could be noticed [START_REF] Bubnicki | Pattern Recognition Algorithms for Simple Markov Chains[END_REF] that using assumption about independence of input observations and that a priori distributions depends on recent a priori distribution and current transition matrix, the expression for dependent risk functional (5) could be calculated using following recursive procedure:

N N N N N P D δ δ ) ( 1 1 1 ⋅ ⋅ = + + + u , 1 1 1 1 ) ( p D ⋅ = u δ (5)
where

T M N N N N ] ... [ ) ( ) 2 ( ) 1 ( δ δ δ δ = and             = ) | ( 0 0 0 ) 2 | ( 0 0 0 ) 1 | ( ) ( M P P P D N N N N N N N N u u u u L M O M M L L . ( 6
)
Hence, to make a decision for a given input in N th time step it is enough to have matrices D N (u N ), P N , and the vector δ N-1 . It is important due to the learning stage.

General methodology

However, from the computational point of view it is important to have as small number of features as possible with smallest loss of information. Therefore, the number of inputs could be decreased to minimum. For example, a student could be described by a mean value of grades, number of courses he attended. Then students could be clustered into groups and each group contains students who have similar descriptions. One group is so called context [START_REF] Chen | Integrating Context-Aware Computing in Decision Support System[END_REF], [START_REF] Palmisano | Using Context to Improve Predictive Modeling of Customers in Personalization Applications[END_REF].

Furthermore, in real situations the probability distributions are unknown. Therefore, a learning stage is needed. During the learning process the distributions are estimated based on training sequence, (u n ,j n ), n = 1,2, ..., N. It is worth to notice, that at each n there could be more than one observation, (u n,k ,j n,k ), k=1,2,...,K. However, in real environments the non-stationarity occurs which means that the estimation cannot be made using all observations, e.g. users' preferences which evolve in time [START_REF] Webb | Machine Learning for User Modeling[END_REF].

Thus, following general methodology could be propose: 1. Divide users descriptions (demands) into clusters (each cluster is called a context) using some chosen clustering algorithm (e.g. k-Means).

Learning Stage:

i.) At the beginning estimate p 1 , D 1 (•) ii.) At each learning step N>1 estimate D N (•), P N , and δ N-1 .

Classification stage:

i.) For given user's demand find an according context (e.g. using 1-Nearest Neighbour classifier). ii.) For given user's context make a decision using (5).

Incremental learning algorithm

As it was mentioned, in many real-life situation the non-stationarity occurs. Therefore, an adaptive estimation methods have to be proposed. Such adaptation to changes could be made by applying incremental learning algorithms [START_REF] Tomczak | Bayesian Classifiers with Incremental Learning for Nonstationary Datastreams[END_REF].

In presented approach a following method of estimation is proposed. Let us assume that at each learning step there are K observations which come in a data stream (sequence). It means that k th observation occurs before (k+1) th observation. Then, the probability distribution is estimated using a frequency matrices (matrices of the same sizes as D N (•), P N , and δ N-1 ) but consisting frequency of of responding input or decision. Let denote such matrix as V D,N , V P,N , V δ,N-1 . Then, the learning algorithm could be as follows 1 :

1. Take a new observation and initiate all matrices.

Update the appropriate frequency matrix by a new

k th observation, V N (k) = a N •V N-1 + V N (k-1)
+ W k where W k is a matrix with ones in proper places for observation's decision and context (determining row and column), and zeros -otherwise, a N ∊ [0,1] is a forgetting factor.

Estimate probability distribution using appropriate frequency matrix. If

there is a new observation, then go to 1.

The key issue is to fix a proper value of the forgetting factor because it affects the estimation. It can have a different value at each learning step or be constant.

Experimental study

Presented approach is checked on the example of the educational platform at WTU. The platform is dedicated to students service and enables e.g. sending applications, signing up for courses, checking past grades, and so on. Therefore, the educational platform could be seen as a SOS in which different services could be distinguished.

However, in this work the proposed approach is used in the enrolment service. At WTU there are different services for enrolment: faculty enrolment, specialisation enrolment, sport enrolment, and foreign languages enrolment. Each of mentioned enrolment needs reading about tens or even hundreds of course descriptions which is a big waste of time both for student and system. In the experiment on the example of a enrolment for foreign language it is to shown how proposed approach could be used.

Experiment details and results

In the experiment the state in the Markov chain is associated with the language and its level: English A1, English A2, English B1, English B2, English C, German A1, German A2, German B1, German B2, German C, nothing. At WTU there are more languages available (e.g. Korean, Japanese, Czech, Italian, Swedish) but for transparency of the experiment it was limited to 11. The nothing means that student is not signed up for any course.

Moreover, it is assumed that each student is described by a following vector: number of all courses, number of exercises, number of laboratories, number of lectures, mean of grades, variance of grades, mean of grades from lectures, mean of grades from exercises, mean of grades from laboratories, number of fails (grade F), number of excellents (grade A). To generate students a real logs from the educational platform was used. This dataset was used also in the previous works [START_REF] Prusiewicz | Services Recommendation in Systems Based on Service Oriented Architecture by Applying Modified ROCK Algorithm[END_REF], [START_REF] Sobecki | Student courses recommendation using Ant Colony Optimization[END_REF].

Then the problem is as follows: Propose a student the most appropriate language and the level if recently she/he was enrolled to a course X.

The experiment was conducted for 2 semester (one semester -one learning stage). In the semester one student signs up for a course after another.

The methodology of the experiment was following: 1. Using the real dataset a context was established using k-Means clustering algorithm. (There are 3 clusters). 2. Student is generated due to the appropriate probability distribution2 and her/his description is mapped with the context id by using 1-NN classifier. 3. Initial state and next states are generated due to (5). 4. Incremental learning algorithm with a forgetting factor a is applied. It was assumed that at the first learning stage (first semester) a 1 = 0. Following methods were compared (for 100 students during a semester, 2000 students, and 5000 students):

• Markov chain (MC) model with Bayesian inference (5) with a 2 ∊ {0, 0.1, 0.25, 0.5}. • Bayesian model (no assumption about Markov chain, Bayes) with Bayesian inference with a 2 ∊ {0, 0.1, 0.25, 0.5}. • Random -list of states is given randomly. Additionally, beside classification accuracy another performance index was used, called position index. This criterion calculates the difference of real state in the sorted due the probabilities list of decisions from the first position. The experiment was conducted on 10 generated datasets (decisions about enrolments) and the results are a mean values. Classification accuracy and position index are calculated as a mean of two semesters. The results are shown in the tables 1 (classification accuracy) and 2 (position index -the lower the value the better).

Discussion

First of all it has to be said that analysing classification accuracy for Markov model and Bayesian model no statistical difference could be noticed. In the situation with 100 students in semester one and two the Markov model performs worst because of its bigger size of matrices. However, in case of 2000 and 5000 students Markov model gave slightly better results. On the other hand, analysing position index Markov model was a little bit better than Bayesian model and the mean value of std shows that Markov model behaves in more stable way (difference of about 0.15). Furthermore, applying Bayesian inference allows to propose around 1/3 of students proper decision or, in average, proper decision at 2 nd or 3 rd position in the list. In comparison to random method it is clear that Bayesian approach outperforms random personalisation. Moreover, the usage of the incremental algorithm gives slight improvement in the quality of classification accuracy and position index.

Concluding, presented experiment is conducted only on two semesters. Probably, carrying out the experiment on 10 semesters should show that applying Markov chain model is more appropriate.

Final remarks

In this paper the general approach using Markov chain model and Bayesian inference for personalisation was proposed. The presented method was evaluated on the partially real-life data take from the educational platform. The problem was formally stated and the general methodology was proposed.

In the future more attention should be paid in developing the whole system with specified, co-working modules. Moreover, given methodology should be applied in the existing educational platform as the additional feature.

Furthermore, the Markov model should be compared with other methods, e.g. using rules-based knowledge representation, or ensemble classifiers. And the higher order Markov chains are supposed to be considered.

Besides, the research ought to be conducted on bigger amount of data, including logs of users containing whole history, e.g. history of all enrolments.

Table 1 .

 1 Classification accuracy (mean value and standard deviation) for compared models.

			Mean			Std
		100	2000	5000	100	2000	5000
	MC 0.0	0,23	0,283	0,285	0,043	0,005	0,004
	MC 0.1	0,229	0,283	0,285	0,042	0,005	0,004
	MC 0.25	0,229	0,283	0,285	0,042	0,005	0,004
	MC 0.5	0,228	0,282	0,285	0,04	0,005	0,004
	Bayes 0.0	0,245	0,281	0,281	0,046	0,009	0,0055
	Bayes 0.1	0,245	0,281	0,281	0,047	0,009	0,0055
	Bayes 0.25	0,244	0,281	0,281	0,047	0,009	0,0055
	Bayes 0.5	0,244	0,281	0,281	0,047	0,009	0,0055
	Random	0,095	0,093	0,090	0,015	0,004	0,001

Table 2 .

 2 Position index (mean value and mean value of standard deviation) for compared models.

			Mean value of position index		Mean value of std of position index
		100	2000	5000	100	2000	5000
	MC 0.0	2,91	2,46	2,42	0,713	0,597	0,583
	MC 0.1	2,85	2,45	2,42	0,662	0,592	0,581
	MC 0.25	2,86	2,45	2,42	0,660	0,592	0,581
	MC 0.5	2,87	2,43	2,42	0,665	0,591	0,581
	Bayes 0.0	2,81	2,53	2,53	0,703	0,745	0,756
	Bayes 0.1	2,78	2,53	2,53	0,636	0,741	0,754
	Bayes 0.25	2,78	2,53	2,53	0,636	0,741	0,754
	Bayes 0.5	2,76	2,53	2,53	0,639	0,740	0,754
	Random	4,95	4,97	6,02	0,230	0,277	0,356

Because the algorithm is the same for all matrices, therefore indexes by V are skipped.

Features number of fails and number of excellents were drawn from discrete distributions and all other -Gaussian.

Acknowledgments. The research presented in this paper has been partially supported by the European Union within the European Regional Development Fund program no. POIG.01.03.01-00-008/08.


João Paulo Barbosa

Paulo Barros

Franklin Ramalho
email: franklin@dsc.ufcg.edu.br

Luís Gomes

Jorge Figueiredo

Filipe Moutinho

André Anikó Costa

Aranha

João Paulo Barros

Anikó Costa

André Aranha

SysVeritas: A Framework for Verifying IOPT Nets and Execution Semantics within Embedded Systems Design

Keywords: Embedded Systems, Petri Nets, Maude, Verification

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

Introduction

It is well accepted that models offer one of the best choices to deal with the development of complex systems [START_REF] Gomes | Behavioral Modeling for Embedded Systems and Technologies[END_REF]. In particular, models improve the communication between developers and customers. However, much more is expected from a single model. For example, in the embedded systems domain, one has to specify the system in an unambiguously way and enable sophisticated system analysis. Due to this fact, and in order to increase consistency, most of the currently accepted model-based development techniques are based on formal models [START_REF] Sgroi | Formal models for embedded system design[END_REF], these models are able to precisely represent the semantics of computation and of concurrency.

Several distinct modeling formalisms, supporting the model-based development attitude [START_REF] Gomes | Formal Methods for Embedded Systems Co design: the FORDESIGN Project[END_REF], have already proved their adequacy for embedded systems design. With Petri nets [START_REF] Claude | Petri Nets for Systems Engineering[END_REF] as a system specification language we get the advantages of its strong mathematical definition and its well defined and precise semantics, enabling the support for simulation, state space generation, and model-checking techniques for verification purposes.

From another point of view, Petri nets are suitable for the definition of automatic model transformations, in order to obtain models at different levels of abstraction. For example, the FORDESIGN project [START_REF] Gomes | Formal Methods for Embedded Systems Co design: the FORDESIGN Project[END_REF] has obtained several interesting results by reusing the benefits from transformations involving models defined using a nonautonomous Petri nets class entitled Input-Output Place/Transition nets (IOPT nets) and the respective code generation for several languages and platforms. Currently, the MDA-Veritas initiative [START_REF] Barbosa | Checking Semantics Equivalence of MDA Transformations in Concurrent Systems[END_REF] has been proposed as one alternative in order to reuse the state of the art in model based development provided by the FORDESIGN project, shifting the focus to MDA (Model-Driven Architecture) [START_REF]OMG: Model-Driven Architecture[END_REF] thus obtaining improvements in the verification of IOPT nets as formal models.

Concerning the aforementioned context, the expected results are highly dependent on the semantics of the IOPT models. Since these models are expected to be built as Platform Independent Models (PIMs), they should be properly transformed into Platform Specific Models (PSMs) for mapping to the corresponding code. However, this automatic generation must take into account the platforms architecture and definition. Thus, an important question arises: how to obtain the verification goals, given the number of existing platforms? Models need to make sense in the corresponding environment they are inserted. Thus, since each specific platform has its own concepts and execution semantics, models should be able to represent these characteristics.

In Section 2, we will discuss about the environmental impacts of deploying a semantic model for PSMs, avoiding the early use of electronic devices. We explore the previously defined semantics for IOPTs, briefly recapitulated in Section 3. In Section 4, an executable formal model upon which deploying environmental characteristics will be represented. Moreover, through a running example in Section 5, the analysis for static properties, synthesis of executable sequential implementations, automated distribution and dynamical behaviors following model checking techniques are presented as one of the main points of this work. The approach is fully supported by MDA standards and tools and takes advantage of the suitability of Maude, and its metalanguage capabilities, for reactive systems modeling and execution. Finally, Section 6 discusses our final view about this gap between the concepts available in a design language, in this case the IOPT nets, and those available in platforms, as well as the analysis and verification tools. We expect that designers should be able to work with the domain-specific abstractions such as signals and events, so the knowledge required to map them into a platform must be provided automatically.

Contributions to Sustainability

We expect that the greater the effort to tackle the question of representing platformspecific characteristics in verification models, the greater will be the obtained benefits. More specifically, the creation of formal models eases the simulation and verification at several levels, hence decreasing the necessity of dealing with hardware devices at early stages, avoiding the use of additional electronic devices that would hardly be recycled, saving energy, increasing the reliability level, and reducing costs. Hence, our approach constitutes a technological innovation that contributes to sustainability. To that end, we have identified four distinct contributions: (i) the use of more reliable development flows, due to the consolidation of a formal MDA approach, (ii) savings in costs due to hardware, energy, specialized engineer views (for the software developer), more abstract specifications, and correctness, (iii) the use of specific executable semantics specifications available in several kinds of platforms and (iv) at the practical level, the facilitation of rapid prototyping for reliable control and domotic systems.

IOPT Nets and MDA-Veritas

Petri nets are a well-known set of formal languages with a common graphical representation, particularly suitable for the visual modeling of concurrent systems. The class of Input Output Place Transition nets extends the class of Place/Transition nets (P/T nets) with non-autonomous constructs allowing the modeling of controllers connected to the environment through signals and events. The respective IOPT semantics that interest us was already presented elsewhere [START_REF] Gomes | Formal Methods for Embedded Systems Co design: the FORDESIGN Project[END_REF]. Next, we briefly present the main characteristics of this semantics.

IOPT nets add several annotations to the P/T nets nodes and net modules. More specifically, transitions can have associated input and output events, as well as a guard that is a function of the input signal values. Additionally, each conflict is resolved through the addition of a priority annotation to each transition. Places have a bound annotation specifying the maximum number of tokens in each place, which can be of major relevance when automatic code generation is considered. They also have a conditional external action on output signals. The respective condition is a function of the place marking. All signals and events are defined at the net module level. Next we briefly present the IOPT nets semantics.

The IOPT nets have maximal step semantics: whenever a transition is enabled, and the associated external condition is true (the input event and the input signal guard are both true), the transition is fired. An IOPT net step is maximum when no additional transition can be fired without generating an effective conflict with some transition in the chosen maximal step. Therefore, we define a IOPT net step occurrence and the respective successor marking. The net evolves through the firing of successive maximal steps. The synchronized paradigm, used in this work, also implies that the net evolution is only possible at specific instants in time named tics. These are defined by an external global clock. An execution step is the period between two tics.

System Modeling

In [START_REF] Barbosa | Semantic equations for formal models in the model-driven architecture[END_REF] we have proposed a formal approach named semantic equations to extract formal models from the syntactic constructors, i.e. the metamodel. There, we have approached Petri nets models involved in model transformations aiming to ensure semantics preserving properties. Here, we are interested in reusing the semantic equations approach but for another purpose: to extract the state space as the semantic representation for the IOPT nets formalism. The understanding of this state space in terms of algebraic structures over graphs has the potential to unify several views of semantics formalisms. This structure can be manipulated by formal methods tools employed in several sorts of analysis and satisfying some requirements that are hard to solve. As an example, several formal tools could handle the state space explosion problem, by providing an on-the-fly way of reasoning over occurrence graph structure of the Petri net model if it is seen as this algebraic structure.

As in [START_REF] Barbosa | Semantic equations for formal models in the model-driven architecture[END_REF], this ordinary graph has a set of nodes as a commutative monoid S+ generated by the set S of places and the empty marking as identity. The sum + of transitions represents the parallel firing of transitions, and the sequential firing of transitions is denoted by the operation. By the guarantee of the closure property, we are representing computations through simple transitions.

The concept of semantic domain and semantic equations gives us a guarantee that the metamodel presented in [START_REF]Ecore Representation for Extending PNML for Input-Output Place-Transition Nets[END_REF] will provide models able to satisfy the desired executability requirement.

From IOPT Models to Algebraic Specifications

.. 6 rl [turnLghtOn]:{DARK}+[SW10]+(waitngDark)=>{idle}+[noSignal]+(light). 7 endm

The following code fragments represent a basic template for the translation of an IOPT model to a rewriting logic specification in Maude. We use regular grammar constructors for specifying the possible number of elements in a model. Line 1 is the declaration of the module able to represent the IOPT net. Line 2 defines the existing sorts for this specification. From line 3 up to 5 we have the declaration of the names for sorts Event, Signal, and PlaceMarking respectively, and the corresponding identity operations, i.e. idle, noSignal and empty. Finally, line 6 is the basic operation of combination of PlaceMarkings that follows the principles of commutative monoid Petri nets representation with the corresponding properties. Thus, the Event, Signal and PlaceMarking elements are composed, in lines 7 up to 9, in EventSet, SignalSet and NetMarking respectively. The composition of these major structures is called a IOPT type at line 11.

Finally, the transitions are represented as rewrite rules from one IOPT configuration (line 12) to another IOPT configuration (line 13).

Maximal Step Semantics

In order to give a semantic representation of Petri nets having the maximal step semantics in a translational way, we need to change the Maude's original interleaving semantics. The main structure of the modules produced in this activity is depicted in Fig. 3. Starting from an IOPT_SYSTEM_CONFIGURATION, our choice relies on the fact that the Maude system contains a predefined module called CONFIGURATION for defining a denotational semantics and a META-LEVEL for redefining the operational semantics. The components from the denotational view were presented in the previous section, through the definition of an IOPT_NET_CONFIGURATION. Here, we focus on the operational view. A the META-LEVEL, we can find operators that represent terms and modules. From this level, the module also supplies efficient descent operations reducing the computations from the metalevel to the object level. Thus, operations such as metaApply, that matches the term with lefthand side of the rule, apply the rule at the top of the term and returns the metarepresentation of the term, can be used in order to take the application of all enabled transitions at a single step. Inheriting from the META-LEVEL, we defined a module META-PETRI-NET able to represent the main structural operations and rules from a Petri net. Finally, the MAXIMAL_STEP redefines the execution semantics of the Maude system for this domain, according to our specification. An excerpt of this module is shown in the following code fragment. Lines 1 and 2 of this excerpt present the basic definition of the maximal-step rule for IOPT nets. A Term t is rewritten of the application of all applicable rules for the given module. Lines 3 and 4 detail the match of a given rule if applicable and the call of the built-in operation metaApply from the Maude META-LEVEL module. Line 5 represents the can in which all applicable rules were already applied, returning the current term.

Example

We have applied this formal specification in different scenarios. In order to illustrate it in practice, we show in Fig. 4 an application under development that supports domotic control systems. For simplification purposes, we have chosen a model with sensors that controls: (i) the arrive and leaving of strange people; (ii) the detection fire and its disabling; (iii) the alert for an unexpected situation and its stop; (iv) the detection of darkness and turning on the lights; and (v) a central enabler that establishes the priorities of these independent actuators.

From this model, we can automatically extract a translational semantic representation able to be employed in the Maude system. The following excerpt illustrates the single transition disFireAlert in Fig. 4 generated automatically from the previously described model. This means that the transition disFireAlert, which represents disabling the alert of fire, has dependencies in its firing from the existence of the event DISABLE and the signal SW4 (representing the action of pressing the switch of number 4 in the board that contains a PIC microcontroller) and having tokens in the places enabled and firing. As result, because of the act of consuming, no event and no signal are available in the system, the matched tokens were removed and new tokens for the places waitingForFire and waitingToEnable were produced. This represents a Platform-Specific Model (PSM) according to the MDA view.

Simulation and Verification

The system can be simulated through common rewrite commands or generate the state space as an excerpt shown in Fig. 5. We developed a plugin to produce graphical visualizations of IOPT states after simulation and verification by using the GraphViz solution [9]. This represents the case where the system just allows the firing of one actuator that disables according to the priorities one enabled module from many (defined as n=1 in the place waitingToEnable) when sensors were fired concurrently. It shows that from the initial state, it is possible, through the maximalstep semantics, the activation of the events ARRIVE, WARM, ENABLE, STOP and DARK. From this state, excepting the luminosity module that is completely independent, only the ALERT event and its corresponding module, sensors and actuators are able to return to the initial state, depending from the ENABLE actuator.

The developers of the domotic model were interested in the verification of properties such as deadlock freeness, ensuring the correct application of the priorities and logical implications that given the firing of a sensor, the corresponding actuator will take the control and after solving will go back to the waiting state. These kind of primitive properties are automatically derived for the Maude LTL model-checker syntax from the initial system marking as the example that follows:

1 eq initial = {noEvent} + [noSignal] + (waitingForPresence waitingForFire 2 waitingForDarkness waitingAlert waitingToEnable waitingToEnable) ... 3 search in DOMOTICS-IOPT-NET : upTerm(initial) =>! Any:Term .

Fig. 5. Excerpt of the generated state space with the maximal step semantics

There, we have a search for a deadlocked state according to the maximal step semantics. The generated verification code starts by translating the initial marking specified in the ecore IOPT model to the lines 1-2. This initial marking specifies that there are no events and signals and the initial tokens. Finally, line 3 represents the excerpt where the search Maude command is performed by translating the initial marking to its representation at the meta-level (upTerm) and try to find any term with no successors (through the command =>!). After successive refinements during the modeling phase, having several improvements, we get a running model with no deadlocks as show the following Maude output.

1 No solution . 2 states: 36 rewrites: 66 in 0ms cpu (0ms real) (~ rewrites/second)

Conclusions

For the embedded system's development, our contribution comes from the fact that safety and economic concerns require high levels of assurance that the designed model will work as expected in a physical environment. In this sense, the move to the physical implementation represents a huge gap of abstraction. Formal models enable us to solve this problem in an elegant fashion. The gap is filled with an artifact that is a representation of the designed model but can also represent logically most of conditions generated by the environment. Therefore we continue in the MDA lifecycle because this technique also represents the conceptual transformations of PIMs to PSMs according to the MDA view.

The present solution still has some limitations concerning the rigorous formalization of some execution semantics regarding the checking of semantics preservation in model transformations, establishing an equivalence relation between the models. More specifically, and although this does not affect the simulation and verification purposes initially established, the state space in the maximal-step semantics case cannot be fully explored for all partial subsets of events generated by the environment. As future work, we intend to extend the solution for more specific platforms, producing a semantic framework that will bring several benefits for the system's designer before producing a device from a chip layout or deploying code in an embedded platform.

Fig. 1

 1 Fig.1shows the main flow for the use of this solution. It is inserted in a major flow of the MDA-Veritas solution. We have two roles: (i) the modeler that starts editing the model until, throughout MDA transformations, generating the Maude code; and (ii) the verifier, that employs the Maude model in concrete syntax until the final analysis of formulae. The formulae can be proved at the object level by using the Maude LTL model-checker through the definition of predicates.

Fig. 1 .

 1 Fig. 1. Modeling and verification process

Fig. 2 .

 2 Fig. 2. Fragment of a domotics model

5

  ops (PlaceMarking_names)* :->PlaceMarking . op empty:-> PlaceMarking . 6 op __:PlaceMarking PlaceMarking ->PlaceMarking [assoc comm id: empty].

(

  12 rl [Rule_name]:{(Event_names)*}+[(Signal_names)*]+((Place_names)*) 13 => {(Event_names)*}+[(Signal_names)*]+((Place_names)*) . )*

Fig. 3 .

 3 Fig. 3. The modules infrastructure

Fig. 4 .

 4 Fig. 4. Simplified domotics model

mod name-IOPT-NET-CONFIGURATION is

sorts Event EventSet Signal SignalSet PlaceMarking NetMarking IOPT .

ops (Event_names) : -> Event . op idle : -> Event .

ops (Signal_names)* : -> Signal . op noSignal : -> Signal .

rl [maximal-step] : T:Term =>

maxStep(T:Term, applicableRules(T:Term, rules, module)) . ...

eq maxStep(T:Term,(rl X:Term => Y:Term [label(Q:Qid)] .) RS:RuleSet) =

maxStep(getTerm(metaApply(module, T:Term, Q:Qid, none, 0)), RS:RuleSet)

eq maxStep(T:Term, none) = T:Term .

Acknowledgment. This work is supported by the cooperation project funded by Portuguese FCT through the project ref. 4.4.1.00-CAPES, and by Brazilian CAPES through the project ref. 236/09.


Adel Ardalan

Ebru Selin Selen

Hesam Dashti

Adel Talaat

Amir Assadi

Design and Applications of Intelligent Systems in Identifying Future Occurrence of Tuberculosis Infection in Population at Risk

Keywords: Mycobacterium tuberculosis, biomarkers and intelligent decision systems, early detection, Vapnik's statistical learning theory

Tuberculosis is a treatable but severe disease caused by Mycobacterium tuberculosis (Mtb). Recent statistics by international health organizations estimate the Mtb exposure to have reached over two billion individuals. Delay in disease diagnosis could be fatal, especially to the population at risk, such as individuals with compromised immune systems. Intelligent decision systems (IDS) provide a promising tool to expedite discovery of biomarkers, and to boost their impact on earlier prediction of the likelihood of the disease onset. A novel IDS (iTB) is designed that integrates results from molecular medicine and systems biology of Mtb infection to estimate model parameters for prediction of the dynamics of the gene networks in Mtb-infected laboratory animals. The mouse model identifies a number of genes whose expressions could be significantly altered during the TB activation. Among them, a much smaller number of the most informative genes for prediction of the onset of TB are selected using a modified version of Empirical Risk Minimization as in Vapnik's statistical learning theory. A hybrid intelligent system is designed to take as input the mRNA abundance at a near genome-size from the individual-to-be-tested, measured 3-4 times. The algorithms determine if that individual is at risk of the onset of the disease based on our current analysis of mRNA data, and to predict the values of the biomarkers for a future period (of up to 60 days for mice; this may differ for humans). An early warning sign allows conducting gene expression analysis during the activation which aims to find key genes that are expressed. With rapid advances in low-cost genome-based diagnosis, this IDS architecture provides a promising platform to advance Personalized Health Care based on sequencing the genome and microarray analysis of samples obtained from individuals at risk. The novelty of the design of iTB lies in the integration of the IDS design principles and the solution of the biological problems hand-inhand, so as to provide an AI framework for biologically better-targeted personalized prevention/treatment for the high-risk groups. The iTB design applies in more generality, and provides the potential for extension of our AIapproach to personalized-medicine to prevent other public health pandemics.

Introduction

Tuberculosis is a severe lung disease which is responsible for increase of 9.4 million [START_REF]Tuberculosis Control: A short update to the[END_REF] cases per year that will results about 2 millions of patients death [START_REF] Volokhov | Mycobacteria Protocols[END_REF]. Infection will be acquired by inhalation of Mycobacterium tuberculosis contaminated air and/or droplets [3]. Even from slight initial invasion of the agent, infection may lead to latent TB or may lead to primary disease [START_REF] Murray | Tuberculosis: The Essentials[END_REF]. The exact time which is taken from the initial infection until the development of disease varies among individuals. This variation can be attributed mainly to the immune status of an individual. Immunocompromised individuals or individuals which their immune systems are suppressed might likely develop primary disease. The World Health Organization WHO [START_REF]Tuberculosis Control: A short update to the[END_REF], has reported the burden of disease under various circumstances, and in particular, the number of bacteria sufficient to infect an individual [3]. After 2008, the numbers of incidence, prevalence and mortality, are estimated to be at least as 9.4 million incidence cases, 11.1 million prevalence cases and 1.3 million deaths [START_REF]Tuberculosis Control: A short update to the[END_REF]. Since an individual can infect 10-15 individuals [3], from a public health viewpoint, early diagnosis and treatment [START_REF] Hopewell | Tuberculosis: The Essentials[END_REF] is crucial to contain the disease. In particular, technologies for early detection and isolation of an infected individual will play a major role in sustainability of the global population health.

Understanding the molecular mechanisms during the invasion of M. tuberculosis provides valuable insights for the analysis of the biological understanding of the course of infection. Besides, the molecular understanding might lead to development of the necessarily better targeted treatment strategies against tuberculosis [START_REF] Triccas | Use of Fluorescence Induction and Sucrose Counter Selection to Identify Mycobacterium tuberculosis Genes Expressed Host Cells[END_REF]. To serve this purpose, [START_REF] Talaat | The Temporal Expression Profile of Mycobacterium tuberculosis Infection in Mice[END_REF] analyzed MTB infected lung samples of immuno-compromised and immune-competent mice. They found the genes that are expressed (or significantly changed) during early invasion through systematic application of the microarray technology [START_REF] Talaat | The Temporal Expression Profile of Mycobacterium tuberculosis Infection in Mice[END_REF]. As a result of this analysis, they reported the differences between expression profiles in three different environments. In an analogous fashion, but in a different context, microarray technology is used to monitor the changes in M. tuberculosis gene expression during the treatment with antituberculous drug isoniazid [START_REF] Wilson | Exploring Drug-induced Alterations in Gene Expression in Mycobacterium tuberculosis by Microarray Hybridization[END_REF]. According to analysis of expressed genes in presence of isoniazid, researchers are more likely to enhance the drug targets. [START_REF] Behr | Comparative Genomics of BCG Vaccines by Whole-Genome DNA Microarray[END_REF] perform microarray analysis in order to understand the differences between genomic structures of M. tuberculosis and M. bovis and other strains of M. bovis which are also compound of the BCG vaccines. Accordingly, this study has aimed to serve the purpose of developing new and more narrowly-aimed vaccines and/or antituberculosis treatment [START_REF] Behr | Comparative Genomics of BCG Vaccines by Whole-Genome DNA Microarray[END_REF]. [START_REF] Fisher | Microarray Analysis of the Mycobacterium tuberculosis Transcriptional Response to the Acidic Conditions Found in Phagosomes[END_REF], draw attention to the function of the acidification during the immune response through using microarray [START_REF] Fisher | Microarray Analysis of the Mycobacterium tuberculosis Transcriptional Response to the Acidic Conditions Found in Phagosomes[END_REF]. As a result of their analysis, they suggest that, acidification might be a signal to induce the gene expression needed by the bacteria to survive against the immune response cells known as phagosomes.

The discussion above is short, but essentially highlights the critically sparse stateof-art knowledge regarding detection, prediction and treatment of individuals at risk, and in fact, almost all categories of individuals infected by this microbe.

In the following, we report on preliminary progress in design of an intelligent system (Section 3) based on our earlier de novo analysis of gene expression timeseries by novel applications of stochastic signal processing, new clustering algorithms, and dynamics of representations of clusters in an appropriate hyperbolic space.

Contribution to Sustainability

One remark encountered often in modern higher education and research is the significance of close collaboration between young investigators in the computational and the engineering sciences with biologists to provide a fruitful framework for the synthesis of diverse concepts and tools. In this way, integration of hardwaresoftware and 'biological knowledgeware' provides prospects of imminent solutions of myriad challenging biological sustainability research problems through collaborative effort. This research is an illustration of this piece of educational wisdom in this regard. The need for solution of hard biological problems has inspired formulation of a number of challenging problems in scientific computation and engineering.

High performance scientific computation plays a critical role in the 21st century research and engineering design optimization. An important case arises in research on challenging problems in Global Public Health, such as sustainability of protective measures against infectious diseases for humans and animals living across all geographic locations, and as much as possible, under all states of living conditions. This research adheres to the above-mentioned objective of sustainability, utilizing the state-of-art in informatics and engineering. Design of intelligent systems have enriched the modern technological societies, and extending its domain to include the entire global community is inevitable for future protection of life and earth and assurance of a sustainable mechanism to provide a healthy society across the globe. Further, among myriad health and disease conditions, there is a serious risk of faster spread of diseases such as TB by more susceptible sectors of individuals at risk due to a compromised immune systems. Sustainability of global health, therefore, must include effective solutions to prevent infection and more likely death by such individuals. On the ethical side, there are limitations to keeping in isolation individuals-at-risk due to the higher risks of eventual severe or fatal sickness. The research on iTB system provides a viable approach to answer the above-mentioned problems based on effective applications of modern engineering and informatics.

iTB: Intelligent TB Dynamical Modeling

We have developed a conceptual framework for dynamical analysis on the grounds of solid mathematical models and empowered by software engineering viewpoint toward information systems development. This viewpoint ensures reusability of the system as a whole or in part for different applications in different disciplines. In other words, a modular design of the system allows us to (1) simplify extending the system capability while maintaining accuracy, (2) rebuild new configurations on demand -e.g. for different clinical applications, and (3) distribution of massive computation among different processing blocks. The entire setting adheres to the engineering principles for deployment within today's High Performance Computation (HPC) infrastructures, like Computational Grids, and most recently, Computation in Clouds. Data-intensive computing era calls for such adaptable architectures to guarantee the applicability of informatics methods for growing population, increasing healthcare demands and personalized (thus, real-time) medical treatment. 1 shows a global architectural view of iTB. The main part is the empirical network construction module which can work iteratively to gain computational access into higher levels of dynamics of massively complex biological systems.

A brief description of modules is provided below. The Imputation Module tries to fill in the missing pieces of the data, which is a common issue in biological timeseries samples. The algorithm utilized here relies on the widely-used Expectation Maximization (EM) algorithm to find the best candidates (i.e. most probable ones) for missing values. We have adapted a generalized version of EM [START_REF] Schneider | Analysis of incomplete climate data: Estimation of mean values and covariance matrices and imputation of missing values[END_REF] to apply to our problem setting. The next Module also performs preprocessing provides a suitable smoothing of the time-series for mathematical analysis, which requires stable local behavior of the time-series (regarding variations and disregarding spiky results). We have employed piecewise cubic Hermite interpolating polynomial (PCHIP) curve fitting approach [START_REF] Fritsch | Monotone Piecewise Cubic Interpolation[END_REF] and re-sampled the regressed curve to approximately simulate the behavior of the system in an appropriately smoother way to accommodate differentiation and other analytic operations.

The next few Modules are designed to find the multi-level structure of the complex networks. They deal with the twisted behavior that results from the connectionism beyond the "build from the simplest blocks" philosophy. In nontechnical terms, connectionist and other learning-theoretic models are constrained by the nature of the domain of generalization, and the balance between the sufficiency of samples versus overtraining. This implies the requirement of elucidating the interrelationship among samples of different levels of significance for estimating future dynamics from the sample point (among other technical hurdles to overcome the biological complexity in predicting the state of the disease from a sparse sample.) Thus, to capture the dynamics of the disease, hierarchical clustering techniques are employed to build a multi-level structural/behavioral model of interactions inside the system. There are advantages to using hierarchical clustering versus non-hierarchical clustering. A comprehensive analysis of different clustering methods and their applicability will appear in due time. Briefly, hierarchical clustering allows the modeler to take into account various forms of analytic singularities, avoiding the artificial assumption that the data samples are uniformly chosen from a single probability distribution function (pdf). Hierarchical clustering, on the other hand, offers the more realistic assumption that different sub-clusters are samples from several pdf that could somewhat differ from the initial pdf. Such diversity of pdf is expected in biology, due to variation and other complex biological phenomena.

The remaining Modules attempt to render the novel concepts of specialarchitecture empirical networks for topological modeling of complex networks. This approach enables us to use a vast spectrum of solid mathematical analysis tools to reveal invaluable measures of correspondence between components of complicated systems (cf. below for an outline.)

Two of the Modules capture dynamic similarities, and record the migration of different classes of genes with different perspectives towards the inclination of the groups' rate and acceleration change in the course of time. The next Module estimates the probability densities in clusters via the exponential family of models. The exponential family has the unique advantage of being quite flexible to accommodate many deviations from Gaussian models, while still are parameterized via a a finite dimensional Riemannian manifold in the Hilbert space of L 2 -integrable functions. The Riemannian structure alluded above is complete and hyperbolic. Hypothesizing the consistent behavior of related genes in the hyperbolic space, we have measured the distance between the clusters of genes as follows, which we mention in the Gaussian case to simplify the presentation.

The individual clusters are regarded as samples from a probability distribution; e.g. for the sake of a concrete illustration, consider two clusters that are regarded as samples from two pdfs that are normal distributions 2 ~( , )

j j j N m µ σ (1)
So for any time frame, we have a topological representation of the system that evolves in time and shows the behavior of the systems in the frame of the hypothesis. To be able to visually inspect the dynamics of the system, we may argue based on the experimental results that projection of the high-dimensional networks obtained into the very first principal components gives us a good representation of the behavioral model.

Mathematical Methods

First, the solution of the preceding estimate requires transforming the ill-posed problem into a regularized well-posed problem [START_REF] Poggio | Regularization algorithms for learning that are equivalent to multilayer networks[END_REF], [START_REF] Girosi | An Equivalence Between Sparse Approximation and Support Vector Machines[END_REF] and [START_REF] Smola | Form Regularization Operators to Support Vector Kernels[END_REF] Thus, it is desired here to have a well-posed problem regarding estimating the measure on the "function-space". There are well-known methods and more modern ongoing research on regularization, and we shall omit the discussion due to lack of space [START_REF] Poggio | Regularization algorithms for learning that are equivalent to multilayer networks[END_REF] and [START_REF] Smola | Form Regularization Operators to Support Vector Kernels[END_REF] To readers familiar with learning theory [START_REF] Vapnik | The Nature of Statistical Learning Theory[END_REF], the latter problem could be regarded as "Learning the Measure" from the sample of trajectory dynamics (we used the data available to us from the TB-infection of mice) through a controlled iterative scheme. Thus, we proceed to cast the latter in Statistical Learning Theory. Accordingly, we need to have a robust estimate for the error in iterative steps of learning to quantify the approximation error for the posterior mentioned above. Robust error estimates require stability in solution of sampling. To have a wellposedness of the inverse problem for the posterior measure will provide desired levels of stability. In turn, such stability may be used as the basis for quantifying the approximation of inverse problems for functions in a finite-dimensional-space setting. This requires an estimate for the distance between the true and approximate posterior distributions, in terms of error estimates for approximation of the underlying forward problem.

Let µ 1 and µ 2 be two normal distributions with means m 1 , m 2 and standard deviations σ 1 and σ 2 . The computation of distance between them is as follows

( ) 2 1 2 2 2 1 2 1 2 1 4 1 2 2 2 1 2 2 1 m m e σ σ σ σ σ σ - - +     -   +   (2)
This metric is initially defined for two measures µ 1 and µ 2 that are absolutely continuous relative to a measure λ as:

1 2 2 1 2 1 2 1 ( , ) 2
H d d d d d d µ µ µ µ λ λ λ       = -           ∫ . ( 3
)
where ௗఓ ௗఒ terms are Radon-Nikodym derivatives respectively. The definition of (3) does not depend on the choice of λ, so (3) does not change if λ is replaced with a different probability measure with respect to which both µ j are absolutely continuous. For two normal distributions (1), the formula for (3) is (2) and readily could be used, even helped by symbolic algebra, to improve accuracy of the iterative 1 2

( , )

H d µ µ

Learning-theoretic calculations. Control of differences in the Hellinger metric (i.e. d 2 H ) leads to control on the differences between expected values of functions and operators (that admit polynomial bounds). Statistical Learning Theory [START_REF] Vapnik | The Nature of Statistical Learning Theory[END_REF] provides the tools to complete the remaining steps in this approach. There are other computational reasons for choice of the Hellinger metric versus other probabilitytheoretic divergences, say from the family of f-divergences such as the Kullback-Leibler, Wasserstein or other metrics. Among them, one could gain useful estimates of bounds more easily using this metric, such as in

1 2 2 2 1 2 1 2 1 2 || || 2( || || || || ) ( , ) H fd fd f d f d d µ µ µ µ µ µ - ≤ - ∫ ∫ ∫ ∫ . (5)
In turn, such bounds point to design of numerical schemes that allow us to solve the inverse problem and gain control in relating estimates arising during perturbations in the domain and range, respectively. Briefly, in the discussion above, let the integer N denote the iteration count, and correspondingly, the estimates

Ψ N (v;w), µ N , such that ௗఓ ௗఓ బ ~݁ିఅ(௩;௪) , ௗఓ ಿ ௗఓ బ ~݁ିఅ ಿ (௩;௪)
. Then bounds on (6) provides a sequence of improving bounds on (3) that demonstrates when ܰ → ∞, the Hellinger metric approaches zero exponentially in N.

| ( ; )

( ; ) | N v w v w Ψ -Ψ . (6)
For Gaussian densities, these bounds are used to prove that the means and covariance operators associated to µ N , and µ, are close in the Hilbert-space (or the Banach space, in more general circumstances) operator-norms. Therefore, in the approach outlined above, we could arrange for the transfer of estimates from the numerical analysis of forward problems into estimates for the solution of the related inverse problem.

Results and Discussion

In the preceding arguments, the space of exponential probability density functions is a Riemannian manifold, and we need a discrete approximation to capture its metric properties within the prescribed error bound. The discrete approximation is typically expected to be high dimensional (thousands or more). In the case of normal distributions, a number of analytic simplifications are available that allow of us to reduce the dimensionality of the metric model. In particular, the Riemannian metric could be approximated by sampling of the distance data, and the model reduction for the sample agrees with the desired approximation to the Riemannian structure. As expected Singular Value Decomposition provides the direct approach, and in the case of animal models of the disease (murine), the results are obtained as follows. Empirically, we have observed that considering the first three principal components retains about 80-85% of the information content of the network (i.e. the ratio of the first 3 eigenvalues to the total sum), while the dynamical separation of different conditions under study are brilliantly visible. We have studied, as an example, the dynamics of gene expression profiles of TBrelated families in a period of 24 hours after infection. Profiles are recorded every 4 hours and a sliding window approach is used to investigate the differences between the "behavioral associations" among genes. As could be observed in the figures, there is a heart-beat-like pattern between clusters representing the conditions. Each point on the graphs represents a cluster of genes in one of the conditions. The observable dynamics resemble a group of particles in a force field which approach each other and then the repulsing forces cause the system to scatter around.

The importance of the above-mentioned results is that the behavioral differences between the two conditions are observable in the very first stages of getting exposed to the invader. The discrimination between the two groups of points (red crosses vs. blue circles) which correspond to the different conditions, could be observed from these graphs as a growth/shrinkage pattern. From a clinical treatment viewpoint, this is of utmost importance as one could be treated before the infection spreads out of control.

At the time of writing this article, we are making progress in translating the above-mentioned visualization-based observations (the differences between the two cases) into a classification scheme within the 3-dimensional principal component space (please see figures). Our method is based on design of new algorithms that minimizes the empirical risk function (analogous to the Vapnik soft-margin SVM empirical risk functional), using the hyperbolic metric in lieu of the Euclidean distance in the theory by Vapnik and others. While the computations are much more challenging, the mathematical results that must guarantee the existence of a minimum for the risk functional and the desired regularity properties are ensured through extensive mathematical work on analysis of functions on hyperbolic spaces, in particular, bounded sequences of approximations to a minimum converge exponentially (hence even faster than the Euclidean metric) and there is a unique limit point for the sequence, hence a unique minimum. The details will be provided in a forthcoming paper.

Fig. 1 .

 1 Fig. 1. Global architectural view of iTBFig.1shows a global architectural view of iTB. The main part is the empirical network construction module which can work iteratively to gain computational access into higher levels of dynamics of massively complex biological systems.A brief description of modules is provided below. The Imputation Module tries to fill in the missing pieces of the data, which is a common issue in biological timeseries samples. The algorithm utilized here relies on the widely-used Expectation Maximization (EM) algorithm to find the best candidates (i.e. most probable ones) for missing values. We have adapted a generalized version of EM[START_REF] Schneider | Analysis of incomplete climate data: Estimation of mean values and covariance matrices and imputation of missing values[END_REF] to apply to our problem setting. The next Module also performs preprocessing provides a suitable smoothing of the time-series for mathematical analysis, which requires stable local behavior of the time-series (regarding variations and disregarding spiky results). We have employed piecewise cubic Hermite interpolating polynomial (PCHIP) curve fitting approach[START_REF] Fritsch | Monotone Piecewise Cubic Interpolation[END_REF] and re-sampled the regressed curve to approximately simulate the behavior of the system in an appropriately smoother way to accommodate differentiation and other analytic operations.

Figures 2 to 7

 7 show a sample dynamics of the TB genes in the level of data. The axes are all relative to a unit-free representation of the cluster distances as measured in the Heilinger distance. Figures are converted to the 3-dimensional projection metric for visualization purposes. The figures are shown in different zoom levels for clarity purposes.

Fig. 2 - 7 .

 27 Fig. 2-7. These samples of the visualization movie frames show instances of patterns in dynamics of the hyperbolic representation of the analyzed sample of TB gene expressions in time. The dynamics is approximated discretely, then projected to the 3-dimensional reduced model from the hyperbolic space. The first three dominant principal modes capture more than 80% of the information contents in the original hyperbolic space. Once the separations of different dynamic patterns are accomplished in the reduced model, clearly the original data will also demonstrate the separation by considering the inverse-images of the separating hyperplanes.

Fig. 4 A

 4 Fig. 4 A sample of projection of the gene expression dynamic pattern in the reduced model.

Fig. 5 Fig. 6

 56 Fig. 5 Another sample of projection of the gene expression dynamic pattern in the reduced model

Fig. 7

 7 Fig. 7 Another sample of projection of the gene expression dynamic pattern in the reduced model.


Ehsan Shahamatnia
email: e.shahamatnia@fct.unl.pt

Ramin Ayanzadeh
email: ayanzadeh@srbiau.ac.ir

Rita A Ribeiro

Saeid Setayeshi
email: setayesh@aut.ac.ir

Adaptive Imitation Scheme for Memetic Algorithms

Keywords: Optimization, Evolutionary computing, Memetic algorithm, Imitation operator, Adaptive imitation

Memetic algorithm, as a hybrid strategy, is an intelligent optimization method in problem solving. These algorithms are similar in nature to genetic algorithms as they follow evolutionary strategies, but they also incorporate a refinement phase during which they learn about the problem and search space. The efficiency of these algorithms depends on the nature and architecture of the imitation operator used. In this paper a novel adaptive memetic algorithm has been developed in which the influence factor of environment on the learning abilities of each individual is set adaptively. This translates into a level of autonomous behavior, after a while that individuals gain some experience. Simulation results on benchmark function proved that this adaptive approach can increase the quality of the results and decrease the computation time simultaneously. The adaptive memetic algorithm proposed in this paper also shows better stability when compared with the classic memetic algorithm.

Introduction and related works

Increasingly, complex systems in different domains raise challenging problems which cannot efficiently be solved with conventional methods. The quest for a solution to these kinds of problems has led researches to use soft computing techniques, by which they can obtain near optimal solutions [START_REF] Eiben | Introduction to evolutionary computing[END_REF][START_REF] Engelbrecht | Fundamentals of computational swarm intelligence[END_REF]. Genetic algorithm (GA) is one of the earliest and most renowned metaheuristics successfully applied on many real world problems [START_REF] Arabi | Cost-effective allocation of watershed management practices using a genetic algorithm[END_REF][START_REF] Donslund | Energy system analysis of utilizing hydrogen as an energy carrier for wind power in the transportation sector in Western Denmark[END_REF][START_REF] Mckinney | Genetic algorithm solution of groundwater management models[END_REF][START_REF] Reca | Genetic algorithms for the design of looped irrigation water distribution networks[END_REF]. Genetic algorithms, like many other metaheuristics, such as particle swarm optimization, explore large areas of search space and locate local minima in early iterations but slack off in trying to find the global optimum [START_REF] Krasnogor | A memetic algorithm with self-adaptive local search: TSP as case study[END_REF]; this behavior is demonstrated in Fig. 1. Another major problem these metaheuristics face is their instability. Due to their stochastic nature they may produce quite different results in different runs of algorithm [START_REF] Eiben | Introduction to evolutionary computing[END_REF]. Among many contributions made to overcome these problems and aiming to improve their performance, another family of metaheuristics, called memetic algorithms, has attracted many researchers. A memetic algorithm, henceforth called MA, is a hybridization of a global optimization technique with a cultural evolutionary stage which is responsible for local refinements [START_REF] Ong | Memetic Computation-Past, Present & Future Frontier[END_REF][START_REF] Ayanzadeh | Determining Optimum Queue Length in Computer Networks by Using Memetic Algorithms[END_REF]. It is reported in the literature that the marriage between global and local search is almost always beneficial [START_REF] Krasnogor | A memetic algorithm with self-adaptive local search: TSP as case study[END_REF][START_REF] Ong | Meta-Lamarckian Learning in memetic algorithms[END_REF][START_REF] Ong | Classification of Adaptive Memetic Algorithms: A Comparative Study[END_REF][START_REF] Liu | An effective PSO-based memetic algorithm for flow shop scheduling[END_REF][START_REF] Moscato | On evolution, search, optimization, genetic algorithms and martial arts: toward memetic algorithms[END_REF][START_REF] Hart | Adaptive Global Optimization With Local Search[END_REF][START_REF] Bersini | Hybridizing genetic algorithms with hill-climbing methods for global optimization: two possible ways[END_REF]. Memetic algorithms are also known in the literature as Hybrid Evolutionary Algorithms, Hybrid GAs, Baldwinian or Lamarckian Evolutionary Algorithms, Cultural Algorithms or Genetic Local Searchers. As biologists put it, a gene is the unit of heredity in a living organism via which the physiological characteristics such as eye and hair color is passed from parents to offsprings. First used in 1976 by Richard Dawkins [START_REF] Dawkins | The Selfish Gene[END_REF], memes are cultural analogues to genes, in that they are responsible for transmitting the behavioral characteristics such as beliefs and traditions [START_REF] Ong | Meta-Lamarckian Learning in memetic algorithms[END_REF]. Biologists believe that genes in a chromosome are (normally) intact during their life span, while psychologists argue that memes are from a more dynamic nature, and they can improve and change under the influence of the environment they are exposed to. Actually the cornerstone of memetic algorithms is that individuals can improve their fitness by means of imitation.

To implement this concept in the computational framework, several approaches have been suggested such as multi-meme MA [START_REF] Ong | Memetic Computation-Past, Present & Future Frontier[END_REF][START_REF] Krasnogor | Multimeme algorithms for the structure prediction[END_REF], coevolving and self generation MA [START_REF] Krasnogor | Toward truly 'memetic' memetic algorithms: discussion and proof of concepts[END_REF] and multi-agent MA [START_REF] Ullah | An agent-based memetic algorithm (AMA) for solving constrained optimization problems[END_REF] , but the mainstream is to conduct a type of local search around each possible solution's neighborhood within a predefined radius [START_REF] Eiben | Introduction to evolutionary computing[END_REF][START_REF] Engelbrecht | Fundamentals of computational swarm intelligence[END_REF][START_REF] Guimares | Optimization of Cost Functions Using Evolutionary Algorithms With Local Learning and Local Search, Magnetics[END_REF]. In the rest of this paper we will consider GA based memetic algorithm but the results can be extended to the most MAs as well. Regarding the interaction of genes and memes in the memetic algorithm evolution process, two main strategies can be witnessed; in memetic algorithms based on Lamarckian theory, the results obtained from imitation operator overwrite the genes. In other words, individuals are strongly affected by the environment and the change in their aptitude is disseminated via the alteration in their genes. In Baldwinian memetic algorithms genes and memes are stored separately [START_REF] Ong | Classification of Adaptive Memetic Algorithms: A Comparative Study[END_REF] as illustrated in figure 2. In this strategy, genes are utilized in the reproduction and memes are utilized in the process of choosing the chromosomes. Hence, among chromosomes with similar fitness, those who have better solutions in their neighborhood have higher chance of reproduction [START_REF] Ong | Classification of Adaptive Memetic Algorithms: A Comparative Study[END_REF]. It is well established that pure genetic algorithms cannot easily fine tune the solutions in a complex search space [START_REF] Krasnogor | A memetic algorithm with self-adaptive local search: TSP as case study[END_REF], whereas if they are coupled with an individual learning procedure capable of performing local refinements the results will be promising [START_REF] Freisleben | New Genetic Local Search Operators for the Traveling Salesman Problem[END_REF][START_REF] Mariano | Arbitrarily large planar etsp instances with known optimal tours[END_REF][START_REF] Ishibuchi | Balance between genetic search and local search in memetic algorithms for multiobjective permutation flow shop scheduling[END_REF]. Benefitting from both Darwinian evolution and domain specific heuristics, memetic algorithms not only balance the generality and problem specificity, but also render a good methodology to trade off between the exploration capacities of underlying GA and the exploitation capabilities of local search heuristic. In the other hand this translates into more computations, more time and loss of diversity within population [START_REF] Araujo | Constrained Optimization Based on Quadratic Approximations in Genetic Algorithms[END_REF]. As an example, with a classic hill climbing algorithm, for a chromosome with n genes if each gene can take 3 different values, then in each iteration 3 n different neighbors must be checked. To ameliorate these drawbacks, we first propose an adaptive imitation scheme and then we study the performance of variations of hill climbing search strategy in conjunction with the adaptive imitation.

Beside their strengths, memetic algorithms also have their limitations. This paper addresses two major issues related to memetic algorithm design and performance. First issue is designing an appropriate local search scheme for the problem at hand. For some problems it is very difficult to define a neighborhood mechanism for a point in the search space. Furthermore, local search strategies usually bear high computational time [START_REF] Gallardo | On the Hybridization of Memetic Algorithms with Branch-and-Bound Techniques[END_REF]. Second issue is setting the imitation operator parameters. In this paper first we propose an adaptive scheme to define the intensity of imitation operator. Then the effect of this adaptive scheme has been studied on different variations of hill climbing local search.

The rest of the paper is organized as follows. In Section 2, we discuss the technological innovation of this paper for sustainability. In section 3, the proposed adaptive imitation operator is evaluated and different local search strategies are studied. The simulation results are also discussed in this section. Finally, Section 4 concludes this paper.

Technological Innovation for Sustainability

This work impacts on sustainability when we are challenged by finding optimal solutions for increasingly complex systems (mostly in terms of dimension of the search space) in different domains, which cannot be efficiently solved with conventional optimization methods. These problems can range from water resource optimization, control and optimization of renewable energy systems, supply chain management and risk management in sustainability proactive strategies. The quest for a solution to these kinds of problems has led researches to use soft computing techniques, by which they can obtain near optimal solutions [START_REF] Arabi | Cost-effective allocation of watershed management practices using a genetic algorithm[END_REF][START_REF] Donslund | Energy system analysis of utilizing hydrogen as an energy carrier for wind power in the transportation sector in Western Denmark[END_REF][START_REF] Mckinney | Genetic algorithm solution of groundwater management models[END_REF][START_REF] Reca | Genetic algorithms for the design of looped irrigation water distribution networks[END_REF]. Memetic algorithms are a class of meta-heuristic algorithms which combine a global optimization technique with a cultural evolutionary stage responsible for local refinements [START_REF] Eiben | Introduction to evolutionary computing[END_REF]. In this paper, we propose a novel adaptive memetic algorithm, where the influence factor of the environment, on the learning abilities of each individual, is set adaptively. This translates into a level of autonomous behavior. In summary, we first propose a memetic algorithm, which includes an adaptive imitation scheme, and then we study the performance of variations of hill climbing search strategy in conjunction with the adaptive imitation. With the improved performance achieved by proposed memetic algorithm we can solve the above mentioned problems more efficiently and achieve better results.

Research contributions and simulation results

As mentioned earlier, the structure of local search algorithms imposes a high computational cost, which is exponential. Simulation results verify that the performance of memetic algorithms is highly dependent to the neighborhood size. For example in a problem with continuous values for genes, smaller neighborhood radius (up to a threshold) leads to better solution, but the price is increased number of evaluations and much computational time. Hence, a dynamic imitation rate in which the neighborhood size is adjusted based on the iteration can guide the searching strategy more wisely. To calculate the imitation rate in continuous search space, equations below are presented:

(1)

a T t rate imitation       - - = Ι 1 1 _ _ (2) ( )         - - = ΙΙ aT t e rate imitation 2 1 2 1 _ _
where t is the current iteration, T is the maximum number of iterations and a is the adaption coefficient. Imitation rate diagram for different values of a is illustrated in figure 3.

The metaphorical concept suggests that imitation can be feckless in some situations, i.e. impulsive imitation misleads the evolution process. Each person during his life span benefits from the interactions with the environment as means to improve his competency in the society, but the scale to which he is influenced from the environment usually decreases as he grows older. In the proposed imitation adaption scheme, diversity of the population is estimated by the iteration number. In the first iterations population should be as diverse as possible to cover all searching space, but achieving final iterations diversity of solutions decreases. In this scheme the imitation rate of individuals decreases along the evolution process. This will not only reduce the unnecessary evaluation of neighbors, i.e. increase the speed, but it also regulates the exploitation rate and leads to improved solutions. Hill climbing local search has been one of the most used canonical imitation operators, suitable for continuous search space [START_REF] Ong | Classification of Adaptive Memetic Algorithms: A Comparative Study[END_REF][START_REF] Liu | An effective PSO-based memetic algorithm for flow shop scheduling[END_REF][START_REF] Bersini | Hybridizing genetic algorithms with hill-climbing methods for global optimization: two possible ways[END_REF]. In the following section different variations of hill climbing algorithm is evaluated and the simulation results are discussed. In these simulations the imitation rate is calculated adaptively using equation 1 and equation 2. To make comparisons we use Ackley test function. Numerous local minima in the search space of this function have made it one of the popular bench mark functions for the assessment of evolutionary algorithms [START_REF] Eiben | Introduction to evolutionary computing[END_REF][START_REF] Guimares | Optimization of Cost Functions Using Evolutionary Algorithms With Local Learning and Local Search, Magnetics[END_REF]. Equation 3 represents the Ackley function in n-dimensional space:

(3)

( ) ( ) π 2 , 2 . 0 , 20 ) ( 1 1 2 cos 1 1 = = = + + ∑ - ∑ - =                 - = = c b a e a e ae X f n i i n i i cx n x n b
Taking into account the stochastic nature of evolutionary algorithms, for each simulation we run the algorithm 30 times and get the best solution and average of the solutions to compare the algorithm performance. Variance of the solutions is also presented to assess the stability of algorithms; smaller variance indicates that solutions obtained from different runs of the algorithm are closer to the average index. Time index indicates the average time a single run of algorithm takes. For all memetic algorithms used in our study, properties of GA part are the same.

Complete hill climbing imitation--In this algorithm all possible neighbors of a chromosome are generated and evaluated and the best neighbor replaces the current chromosome. The search continues until there is no any neighbor better than current state of the chromosome. The simulation results obtained from 30 runs of the algorithm are provided in table 1. To calculate the imitation rate equation 1 and equation 2 are used. Adaption coefficient is set empirically for each problem space. Since the number of iterations required for imitation operator is not predictable and the computational cost is exponential, running time for this algorithm differs in various runs and even the algorithm may end up acting like a complete search. 2First-best hill climbing imitation--Contrary to the complete hill climbing search, this algorithm stops the search procedure as soon as a neighbor better than current chromosome is found. In the worst case first-best hill climbing search is like complete hill climbing search, but the simulation results in table 2 show that the running time of this algorithms is improved compared to the previous one. Single-step hill climbing imitation--This imitation operator considers a hill climbing search with only a single step (iteration). After all the neighbors of a chromosome are evaluated, best neighbor (if any) replaces the current chromosome and search stops. Since each imitation operation is a single step forward, the computational time of the algorithm in different runs is the similar. Simulation results are provided in table 3. Table 4 provides the simulation results to compare an enhanced GA with elitism, standard GA based memetic algorithm (MGA) and single-step memetic algorithm (SSMGA) together. It can be seen that the memetic algorithm with single-step hill climbing search and using the adaptive imitation outperforms other algorithms both in computational time and in the performance.

∞ ---------- ---------- ---------- SSMGA 2

Conclusions and Future Works

Memetic algorithms, being the most famous hybridization method, successfully improve the performance of their underlying global search algorithm. They effectively improve the stability and fine tune the algorithm to converge to the global optimum. Two major issues that MAs suffer from are increased computational time and designing appropriate imitation operator which defines the neighborhood function. Points that are locally optimal with respect to one neighborhood structure may not be with respect to another [START_REF] Smith | Coevolving Memetic Algorithms: A Review and Progress Report[END_REF]. In this paper first we introduced an adaptive scheme to determine the imitation rate. In this scheme the capacity of individuals to be influenced by the environment is reduced as they grow older. By studying different imitation structures, we introduce a single-step hill climbing search strategy coupled with adaptive imitation. Improvements in quality of results and computation time are usually conflicting objectives. Simulation results show that the proposed approach not only improves the efficiency of the memetic algorithm, by fewer neighbor evaluations and less computing time, but it also improves the performance by achieving better results. This verifies that excessive imitation can potentially distract the algorithm from global convergence.

Proposed contributions have been developed for problems with continuous search space, using genetic algorithm and variation of hill climbing algorithm, but the concepts can be extended to discrete domain and using other metaheuristics. Future research will further address the automatic fine tuning of the adaption coefficient.

Fig. 1 .

 1 Fig. 1. Sample performance of a genetic algorithm in minimizing a test function.

Fig. 2 .

 2 Fig. 2. Sample chromosome in Baldwinian MA

Fig. 3 .

 3 Fig. 3. Imitation rate for different values of a. Left: Equation 1, Right: Equation 2

Table . 1

 . . Simulation results for complete hill climbing with adaptive imitation

	Adaption equation	a	Time (s)	Solution average	Solution variance	Best solution
	Eq. 1	3 16.02	4.368 × 10 -8	3.174 × 10 -14	0.000 × 10 -3
	Eq. 2	5 16.02	8.839 × 10 -8	4.329 × 10 -14	3.553 × 10 -15

Table . 2

 . . Simulation results for first-best hill climbing with adaptive imitation

	Adaption equation	a	Time (s)	Solution average	Solution variance	Best solution
	Eq. 1	3 10.23	5.310 × 10 -9	6.672 × 10 -14	0.000 × 10 -3
	Eq. 2	4 10.22	3.199 × 10 -9	5.648 × 10 -14	3.553 × 10 -15

Table . 3

 . . Simulation results for single-step hill climbing with adaptive imitation

	Adaption equation	a	Time (s)	Solution average	Solution variance	Best solution
	Eq. 1	5 6.43	3.165 × 10 -9	4.677 × 10 -17	3.553 × 10 -15
	Eq. 2	2 6.43	1.701 × 10 -11	3.032 × 10 -22	0.000 × 10 -3

Table . 4

 . . Simulation results comparing GA and variations of MA with adaptive imitation

	Algorithm	a	Time (s)	Solution average	Solution variance	Best solution
	GA	-	4.32	4.09 × 10 -6	2.70 × 10 -10	5.58 × 10 -15
	MGA	10 -1	11.2	1.10 × 10 -6	3.61 × 10 -11	5.55 × 10 -15
	MGA	10 -3	285.2	7.49 × 10 -7	9.33 × 10 -12	3.55 × 10 -15
	MGA	10 -6	3892	1.81 × 10 -8	1.67 × 10 -13	2.51 × 10 -15
	MGA	10 -9