Teresa Podsiadły-Marczykowska

Rafał Zawiślak
email: rafal.zawislak@p.lodz.pl

Experimental Verification of the Effectiveness of Mammography Testing Description's Standardization

Keywords: mammogram interpretation, biomedical engineering, ontology, breast scanning, cancer diagnostics

The article presents assumptions, and results of a test of experimental verification of a hypothesis stating that the use of the MammoEdit -a tool that uses its own ontology and the embedded knowledge of mammographyincrease the diagnostic accuracy as well as the reproducibility of mammographic interpretation. The graphical user interface of the editor was similarly assessed, as well as the rules for visualization which assists the radiologist in the interpretation of the lesions' character.

Introduction

Mammography is well founded, effective imaging technique that can detect breast cancer in early and even pre-invasive stage, widely used in screening programs. However, the diagnostic value of mammography is limited by significant and up to 25% high, rate of missed breast cancers. Mammography is commonly seen as the most difficult imaging modality. Computer-aided Detection (CADe) systems aiming to reduce detection errors are now commercially available and are gaining increasing practical importance, but no serious attempts have been made to apply Computeraided diagnosis (CADx) systems for lesion diagnosis in practical clinical situations. The paper presents an assessment of diagnostic accuracy of MammoEdit -an ontology-based editor supporting mammograms description and interpretation. In our opinion MammoEdit fills the gap among the existing solutions to the problem of supporting the interpretation of mammograms. While working on the MammoEdit project, they used the ontology of mammography, created specifically for this need. It was used as a partial set of project requirements for the user interface and database that stores patients descriptions. This role of the ontology of mammography in the MammoEdit editor project was in line with literature indications for the use of ontology in IT systems [START_REF] Podsiadły ; Marczykowska | The role of domain ontology in the design of editor for mammography reports[END_REF].

The Ontology of Mammography

The mammographic ontology was created using Protege-2000 editor and the OWL language. Presented version of ontology includes: comprehensive, standard description patterns of basic lesion in mammography (masses and microcalcification clusters) enhanced with subtle features of malignancy, classes presenting descriptions of real lesions carried out on the basis of completed patterns and grades modeling diagnostic categories of the BI-RADS system.

Radiologists usually agree in their assessments of the diagnostic categories of BI-RADS (Kappa=0.73), when they have to deal with lesions presenting features typical for radiological image of breast cancer, clusters of irregular, linear or pleomorphic microcalcifications, or with spicular massess. The agreement in assessing is significantly lower when it comes to lesions from changes of the 4 BI-RADS1 category (Kappa=0.28). Why does it happen? Admittedly, the BI-RADS system controlled vocabulary contains terms important with specific diagnostic value, but it lacks complete lesions definitions and does not stress the importance of early and subtle signs of cancer. Moreover BI-RADS system recommendations for estimating lesion diagnostic category are descriptive, incomplete and imprecise, they are expressed using different generality levels. The recommendations refer to the knowledge and experience of the radiologist and to typical images of the mammographic lesions. As a consequence, there is a large margin of freedom left for individual interpretation. Those conclusions are confirmed by large study based on 36 000 mammograms [START_REF] Miglioretti | Radiologist Characteristics Associated With Interpretive Performance of Diagnostic Mammography[END_REF], showing wide variability in mammograms interpretation. While creating ontology, ambiguity in BI-RADS system recommendations caused by the lack of data, necessary to express restrictions on allowed values of the ontology classes modeling diagnostic categories of the BI-RADS system. The lack of these classes in the domain model clearly prevents the use of the ontology's classification mechanism to the diagnostic categories' assessment of the real lesions. The important part of this presentation is filling of the gap in imprecisely formulated BI-RADS system recommendations. This work was necessary to define through ontology classes modeling BI-RADS categories.

MammoEdita Tool Supporting Description and Interpretation of the Changes in Mammography

Analyzing reasons for errors it can be generally concluded that potential source of mammograms' interpretation mistakes is uneven level of knowledge and diagnostic abilities of radiologists and their subjectivity. Considering the above, it should be accepted that reduction of the interpretation errors can be obtained by describing mammograms using standardized protocols which include complete, standardized definitions of the lesions taking into account initial, subtle signs of malignancy and indication of feature values qualifying the nature of the changes during the description (potentially malignant or benign). Suggested, by the authors of the article, method of errors' reduction is to support radiologists with specialized IT tool editor -MammoEdit. In the case of the Mam-moEdit, created ontology of mammography was used as a set of partial project needs for the user's interface and for database storing patients descriptions. The main project assumption of the interface was that the description of the mammogram should be done in the graphical mode, using such components as menu and button, while the text of the final description should be generated automatically. The choice of the graphical mode seemed to be obvious, because radiologists are used to images. Pictographic way of recording information enables both fast entry and immediate visualization. Over 300 pictograms illustrating the most important features of the changes in a vivid and explicit way have been created for the need of the editor. The collection of graphical objects (projects of buttons, icons, pictograms etc.) was gathered after the ontological model had been created and it reflects its content. The system of mutual dependence of entering/representing data is also based on the ontology. It let to obtain legibility and explicitness of the interface and blocked entering conflict data (exclusive). Detailed description of the editor can be found in [START_REF] Podsiadły ; Marczykowska | The role of domain ontology in the design of editor for mammography reports[END_REF].

4

The evaluation of the mammogram's interpretation supporting tool

The basic aim of this testing was experimental verification of the hypothesis which states that: "the mammogram's interpretation, using the editor that systemizes the description process and indicates the values of the features that are significant for the accurate interpretation of their nature to the radiologist, increases diagnostic efficiency of the radiologist".

Receiver Operating Characteristic analysis (ROC) has been used to assess testing. Quantitative assessment of diagnostic test performance involves taking into account the variability of the cut-off criteria in the whole variation of the parameter. ROC curves, which are created as a result of drawing the correlation of test's sensitiveness (Y axis) in Specificity function (X axis), are served for this purpose. Sensitiveness (SE) is a test ability to detect disease among sick patients, specificity (SP) is a test ability to exclude disease for healthy patients. The whole field contained under the ROC curve (AUC -Area Under Curve) is interpreted as a measure of test's effectiveness [START_REF] Swets | ROC analysis applied to the evaluation of medical imaging tests[END_REF][START_REF] Zhou | Statistical Methods in Diagnostic Medicine[END_REF]. It has been decided to use Multireader Receiver Operating Characteristic analysis for testing MammoEdit. The mathematic model of the method takes into account typical for radiology sources of variability in image testing assessment and possible correlations among doctor's assessments within the same technique and also possible correlations among assessed techniques. Comparing the effectiveness of different image testing interpretation techniques involves carrying out an experiment, where trial testing consists of: Assessing techniquesimaging, interpretation or medical image processing algorithms; Assessing cases (the collection must include the study of pathological states and negatives); Radiologists who interpret testing cases using testing techniques. These kinds of experiments, marked briefly as "technique x radiologist x case" are based on, so-called, factorial design [START_REF] Obuchowski | Multireader, Multicase Receiver Operating Characteristic Analysis: An Empirical Comparison of Five Methods[END_REF]. The choice of the scale, which was used for assessing test results reflects substantive criteria of the medical problem that is being considered. [START_REF] Zhou | Statistical Methods in Diagnostic Medicine[END_REF]. Diagnostic category of the BI-RADS system matched by radiologists was chosen to represent the assessment of test results. The test which was performed included 80 mammograms that were difficult to interpret (20 negatives and 60 pathologies) The mammograms were described by three different groups of radiologists representing different levels of competence (a trainee, a specialist, an expert). They used two interpretation methods with and without MammoEdit editor. The mammograms came from a public data base DDSM2 (Digital Database for Screening Mammography). Taking under consideration different origin of testing, doctors, and widely varied, uneven level of professional competence, model: random testing/random viewer was used.

The conditions of mammograms interpretation mimicked clinical environment as closely as possible. The mammograms were presented to the radiologists on the medical screens with the MammoViewer application, which function for testing had been reduced to an advanced medical browser, in a dark room. The assessment of mammograms was held during many sessions that lasted from 1 to 3 hours. The speed and time of sessions were under doctors control, the breaks or suspensions were taken on their requests. During one session the doctor performed only one type of mammogram's assessmentwith or without the assisting tool. Changes were assessed on the basis of breast images in two basic projections, without any additional diagnostic projections, clinical trials, patient's medical records or previous images to compare. Before the first mammograms reading, the meaning of scales was explained and, during testing, the scales description was available at the workplace. The minimum time interval between the interpretation of the same testing with and without MammoEdit was estimated on the basis of literature. It is assumed that it should be minimum 4-5 weeks. In the experiment it ranged from 7 to 12 weeks.

The DBM method of experiments assessment MRMC ROC [START_REF] Dorfman | Monte Carlo validation of a multireader method for receiver operating characteristic discrete rating data: factorial experimental design[END_REF] was used to calculate indicators of diagnostics performance, while the PropRoc method [START_REF] Pesce | Reliable and computationally efficient maximum-likelihood estimation of "proper" binormal ROC curves[END_REF], which provides reliable indicators estimation for the analysis of the small trial testing and discrete assessment scale, was used to estimate ROC curves indicators.

The use of MammoEdit to support mammograms interpretation has raised the average diagnostic performance (statistically significant increase, significance level less than 0.05). Detailed results are presented in Table 1. The rise of indicators mentioned above was accompanied by the reduction of the deviation/variation of the standard field under the ROC curve. The biggest increase of the diagnostic effectiveness was observed in case of a specialist radiologist, the smallestin case of an expert radiologist. It proves bad influence of routine on the quality of mammograms assessment. The influence of Mam-moEdit on the variation of radiologists diagnostic opinions was assessed comparing the value of Kappa statistics achieved in the I and II testing period. The variability of opinions was assessed for the pair of radiologists and for the group of all radiologists, separately for microcalsifications, tumors and for both types of pathologies.

The use of MammoEdit for the mammograms interpretations has increased diagnostics opinions consistency for the group of radiologists for both types of changes masses and microcalcification clusters and for each of the pathologies; the biggest increase in consistency has been observed for microcalsifications clusters. MammoEdit has also influenced the increase of opinions consistency in all pairs of viewers and in all groups of changes.

Conclusion

The use of MammoEdit for mammograms interpretation support has raised the average diagnostic performance for the group of radiologists, other ROC analysis indi-cators have also increased. It brings the conclusion, that initially, in the trial testing, the hypothesis of MammoEdit usefulness as a tool for mammograms interpretation support has been confirmed. The level of diagnostic knowledge representation in the ontology of mammography has also been estimated as a satisfactory one. Thorough analysis of radiologists performance, in distinction of the types changes, has proved a difference in supporting tumors and microcalsifications interpretation effectiveness (to the detriment of tumors).

The next result of the experiment is the initial assessment of the logical correctness of the class construction method. It consisted of comparing BI-RADS category of selected testing cases in DDSM data base and the assessment of the grades of the same changes performed using the MammoEdit and automatic evaluation of the changes grade on the basis of ontology classification. The experiment has proved that BI-RADS category assessment of the changes obtained with the help of the ontology's inference module is adequate to the ones in DDSM data base. What is interesting, the use of ME editor has improved asses compliance with the standard, even for the group of expert radiologists. MammoEdit supports just a part of radiologists work, we need to gain some knowledge about mathematical descriptors of all basic types for feature values of the mammograms changes, tumors and microcalsification clusters to apply it in clinical practice and integrate it with the application CAD. We also need to connect changes detection and their verified description, and finally widen the application with the cooperative function with the ontology of mammography in order to receive an automatic verification of the BI-RADS category of the change. The complement of the integrated application functionality is: administrative features (Radiology Information System), data management (data base functions) and security, finally integration with a workstation and PACS (Picture Archiving and Communication System).

Fig. 1 .

 1 Fig. 1. MammoEdit's influence on the diagnostic effectiveness of the radiologists

Table 1 .

 1 The result of MammoEdit editor's influence on mammograms interpretation perfor-

		mance with 3 groups of radiologists
	lesion type	interpretation methods	AUC	SL	AC%	SE%	SP%
	masses	without MammoEdit with MammoEdit	0,610 0,853	0,016	55,0 79,3	20,0 63,3	94,2 97,1
	microcalcifications clusters	without MammoEdit with MammoEdit	0,677 0,730	0,728	71,3 75,9	36,1 44,4	96,1 98,0

These are suspicious changes with ambiguous configuration of feature values, without evidence of malignancy

The material comes from four medical centers in the USA (Massachusetts General Hospital, University of South Florida, Sandia National Laboratories and Washington University School of Medicine).


Mu Lin
email: mlin@math.tsinghua.edu.cn

Jinqiao Yang

Bin Xu

Modeling and dynamic analysis on animals' repeated learning process

Keywords: dynamical modeling, animals' learning, animal training, piecewisesmooth system

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

Introduction

Learning is always a hot issue in pedagogy and psychology research. J. Dewey considered learning as the process of individual activities, from the perspective of Pragmatism [START_REF] Dewey | Democracy and Education: An Introduction to the Philosophy of Education[END_REF]. He proposed "learning by doing", and emphasized the effect of specific situation on the individual learning behavior. However, from the perspective of Constructivism, J. Piaget considered learning as a kind of interaction process between individual and the environment surrounding him, which makes his cognitive structure develop constantly and gradually [START_REF] Piaget | The Psychology of Intelligence[END_REF].

We can see different learning behaviors everywhere. One common adopted way is repeated learning, such like reciting foreign language words repeatedly to remember them, and doing exercises frequently to learn driving. [START_REF] Schmidt | A schema theory of discrete motor skill learning[END_REF]. This process can make some adaptable changes on the mental level. And one of the most important changes is, learners can abstract an "schema" from the operating process of the skill, and this is so called "transfer of learning" [START_REF] James | Theories of the transfer of learning[END_REF].

In modern psychology, learning is such a kind of behavior which is not only owned by human beings, but also by animals. Lots of researchers studied animals' learning behavior. In experiments, E. Thorndike observed the escaping behavior of a cat which was trapped in a box. He argued that the cat keep running aimlessly at first but gradually drop the ways that cannot open the box and finally acquire the right way to escape [START_REF] Thorndike | Animal Intelligence[END_REF]. B. Skinner observed the bar-pressing behavior of a mouse in the box, which can serve itself food, and argued that the animals reacted to the reinforcement in particular ways. And this "reinforcement -reaction" connection is learned by the animals themselves [START_REF] Fester | Schedules of Reinforcement[END_REF].

There are a lot of studies of learning behavior of animals, although the learning contents maybe different, but the way of learning is always repeated learning, for it is widely believed that animals' learning ability is poor. Furthermore, lots of studies seem to tell us that based on this level of learning ability, the animals can perform the transfer of learning like human being. For example, in the Pavlovian classical experiments, the dog acquired the conditioned reflex of "bell ringing -salivating" [START_REF] Chance | Learning and Behavior: Active Learning Edition[END_REF]. So, if we make an animal learn a kind of behavior repeatedly, it indeed can influence the generation of other related behaviors. But if this influence exists, how it works? And how strong it can be? Both of these questions are important because they can help us know the nature of learning deeper. The present study will answer these questions, and try to provide a mathematical model of relationship between these two behaviors.

Mathematical modeling is in psychology and society has rapid development these years [START_REF] Burghes | Modelling with Differential Equations[END_REF]. Social mechanisms are analyzed by analytical and numerical approaches. For example, He et al. found a cellular automata model to describe the learning process [START_REF] He Mingfeng | A cellular automata model for a learning process[END_REF]. Dynamic analysis is used for modeling the problems that have definite dynamic relation or feedback relationship. It has been widely used to solve ecology and demography problems. There are obviously two subprocesses in this repeated learning problem, i.e., the learning are forgetting states. It belongs to the piecewise-smooth problems [START_REF] Kunze | Non-smooth Dynamical Systems[END_REF].

The remainder of this paper is as follows. Section 2 is the mathematical modeling process, and three types of switching rules are enumerated in Section 3. More precise dynamic analysis including some numerical simulation is done in Section 4. Section 5 is the conclusion.

Modeling Framework

We establish a dynamic model in this section first. As mentioned above, for artificial animals' repeated learning process, there are obviously two discrete states, the training condition and forgetting condition, noted as Q = 1 and Q = -1 separately. Note x as the main memory ability, and y as the relevance memory ability, 0 ≤ x, y ≤ 1.

In the training condition, following the logistic model, the dynamic function of the main memory can be described as

ẋ = α 1 x(1 -x).
(1)

While in the forgetting condition, by Ebbinghaus forgetting curve, the dynamic function of main memory ability can be described as

ẋ = -β 1 x. ( 2
)
The variation of the relevance memory ability is related with the main memory ability. In the training condition, it can be described as

ẏ = α 2 (x -l 1 )y(1 -y), ( 3
)
where α 2 > 0. There is one extra term (x -l 1 ) describing the influence of the main memory ability to the related memory ability. More concretely, there is a separation x = l 1 , 0 < l 1 < 1, such that the relevance memory ability increases when the main memory ability x > l 1 and reduces when x < l 1 . Similarly, the equation is

ẏ = -β 2 (l 2 -x)y(1 -y) (4)
in the forgetting condition, where 0

< l 1 < l 2 < 1, β 2 > 0.
The explicit solution of the simultaneous equations ( 1)-( 4) is

x(t) = 1 - 1 1+e α 1 t+c 1 y(t) = 1 - 1 1+(1+e α 1 t+c 1 ) α 2 α 1 e -α 2 l 1 t+c 2 when Q = 1; x(t) = c 3 e -β1t y(t) = 1 - 1 1+e -β 2 l 2 t-c 3 β 2 β 1 e -β 1 t +c 4 when Q = -1, ( 5
)
where c i , i = 1, 2, 3, 4 are arbitrary constants.

In order to get a complete deterministic system, a switching rule between the training and forgetting conditions is needed. Actually, there are four following types, as listed in the following section.

3 Dynamic analysis under different switching rules

Quantitative control switching rule

The quantitative control switching rule contains two preestablished bounds w 1 and w 2 . During the training condition(Q = 1), the system switches to the forgetting condition(Q = -1) when x reaches w 2 . Vice versa. So it makes a periodic learning cycle consist of one training condition and one forgetting condition, See equation [START_REF] Fester | Schedules of Reinforcement[END_REF].

Q(t + ) = -Q(t) if Q(t) = -1 and x(t) ≤ w 1 , or Q(t) = 1 and x(t) ≥ w 2 Q(t) else , (6)
Composed the above switching rule and equation ( 1)-( 4), a deterministic piecewise-smooth system is obtained. Usually, the analysis of piecewise-smooth systems is quite complicated, while our model is a providential one for it has an explicit solution. Assume the point (w 1 , y) changes to (w 1 , f (y)) in a whole training-forgetting period. By computation ,we get f (y) =

Ay (A-1)y+1 , where A = [( w1 w2 ) l2 e w2-w1 ] β [( w1 w2 ) l1 ( 1-w1 1-w2 ) 1-l1 ] α > 0, and α = α2 α1 , β = β2 β1 . So f (y) = A [(A-1)y+1] 2 .
Consider the equilibrium point of map f , i.e. the point y satisfies f (y) = y. The equilibrium point of f indicates the period solution of the whole system, meaning that the collateral memory ability changes periodically with the main memory ability. Solving the equation f (y) = y, we get A = 1 or y = 0, 1. For y = 0 and y = 1 are two trivial solutions, we focus on the case A = 1, then

α(1-l 1 ) ln(1-w 1 )+(αl 1 +βl 2 ) ln w 1 -βw 1 = α(1-l 1 )ln(1-w 2 )+(αl 1 +βl 2 ) ln w 2 -βw 2 Note g(w) = α(1 -l 1 ) ln(1 -w) + (αl 1 + βl 2 ) ln w -βw, then g(0) = g(1) = -∞. The solution of g (w) = 0 is w 0 = l 2 - √ (α+β-βl2) 2 +4αβ(l2-l1)-(α+β-βl2) 2β
. A group of parameters α i , β i , l i , i = 1, 2 guaranteeing A = 1 may exist for w 1 , w 2 satisfied 0 < w 1 < w 0 < w 2 < 1 can be easily chosen. That means f (y) ≡ y, ∀0 < y < 1, i.e., the collateral memory ability can keep a stable circle at any ability.

Otherwise if A > 1, then f (0) = A > 1, so y = 0 is an unstable equilibrium point, while f (1) = 1/A < 1, so y = 1 is a stable equilibrium point. That is to say, under this quantitative control switching rule, orbits start from every initial point except x = 0 will converge to y = 1 after sufficient times of trainingforgetting cycling, i.e. the collateral memory ability will converge to the upper bound. Oppositely, if A < 1, y = 0 is a stable equilibrium point while y = 1 is an unstable equilibrium point. So the collateral memory vanishes no matter how the initial ability was.

Timing control switching rule

Another common switching rule is based on timing control. In this case, each time length of training and forgetting condition is fixed, noted as τ 1 and τ 2 separately. In most cases τ 2 > τ 1 > 0.

The Poincaré map of the system can also be computed explicitly in this case. The system starts from (x 0 , y 0 ) and changes to

x(τ 1 +τ 2 ) = e -β1τ2 1 + ( 1 x0 -1)e -α1τ1 = x 0 • 1 e β1τ2 [x 0 (1 -e -α1τ1 ) -(e -β1τ2 -e -α1τ1 )] + 1 and y(τ 1 + τ 2 ) = M M +N ( 1 y 0 -1) , where M = e -β2l2τ2-β 2 β 1 x1e -β 1 τ 2 , N = (x0) α 2 α 1 e β 2 β 1 x 1 (1+( 1 x 0 -1)e -α 1 τ 1 ) α 2 α 1 e α 2 (1-l 1 )τ 1
, after a period.

So the equilibrium point is x = 0 or x = e -β 1 τ 2 -e -α 1 τ 1

1-e -α 1 τ 1

. For the model located in x ≥ 0, there is only one equilibrium point x = 0 and it is stable, when α 1 τ 1 ≤ β 1 τ 2 . There are one unstable equilibrium point x = 0 and one stable equilibrium point x(0) = e -β 1 τ 2 -e -α 1 τ 1

1-e -α 1 τ 1 . In this case, y(0) = y(τ 1 + τ 2 ), ∀y(0) if and only if M -N = 0. Otherwise the system has two equilibrium points y = 0 and y = 1. If M > N , 0 is unstable and 1 is stable, while if M < N , 0 is stable and 1 is unstable.

Semi-timing control switching rule

Composing the above two rules, the semi-timing control switching rule is also widely used in actual. This rule fixed simultaneously the upper bound (or lower bound) and the forgetting time length (or the training time length). The corresponding Poincare map is:

(1)The system starts from (w 1 , y 0 ), goes through a τ 1 time training condition and reaches (x 1 , y 1 ), then goes back to (w 1 , y 2 ) through a forgetting condition. By computation, we get

y 2 = M M +N ( 1 y 0 -1) , where M = (w 1 + (1 - w 1 )e -α1τ1 ) ( β 2 l 2 β 1 + α 2 α 1
) e α2(1-l1)τ1 , N = e w2 -1)e -α1τ1 )

α 2 α 1 ( 1-w2 e β 1 τ 2 -w2 ) -α 2 α 1 (1-l1) , N = x - α 2 α 1 0 e β2l2τ2+ β 1 β 2
w2(e -β 1 τ 2 -1) . In both two subcases, the dynamic property of y is similar with the one in the timing control switching rule. [START_REF] James | Theories of the transfer of learning[END_REF] The ultimate stable states and bifurcation analysis Although infinite multiple stationary states are existence in some cases of the above analysis, most of them are unstable in fact, for the observation is discrete. More specifically, throughout the above analysis, we assume that the system switches instantaneously whenever the main memorial ability reaches the preestablished value in Section 3.1 and 3.3. Actually, this instant transformation cannot be realized for the estimation of memorial ability must be discrete, furthermore, it costs time.

To describe the influence of the discrete observation, the switching rule can be rewritten as: (take the quantitative control switching rule in Section 3.1 for instance)

Q(t + ) = Q(t) if t n ≤ t < t n+1 -Q(t) else , ( 7
)
where t i , i = 0, 1, ... is a group of switching time series. Usually they are determined mainly depending on a presupposed admissible error ε > 0. And the system switches when the observed variable enters a presupposed neighborhood of the border line. That is to say, if ∆t is the observed interval time length,

t i = k i ∆t, s.t.|x i (k∆t) -w 1 | < ε, or |x i (k∆t) -w 2 | < ε, k i ∈ N.
Under this switching rule, the ultimate stable stationary states can be drawn by numerical simulation, as in Figure 1. Unlike the results given in Section 3.1, there are only finite stable stationary states, means only these states can exist under the discrete observation.

Conclusion

Mathematical model is established to discuss the relationship between the main memory ability and the relevance memory ability. Different switching rules are introduced and under each rule, equilibria and their stability are discussed. Although with some particular parameters, every orbit can be a periodic invariable one, but there can only be finite periodic orbits, due to the discrete observation, as explained concretely in Section 4. The results of dynamic analysis can be used for animal training. Further observational study can be done to demonstrate the models.

Acknowledgements

This study was funded by the National Natural Science Foundation of China(11072274).

- β 2 β 1 1 w 1 + 0 - 1 )

 1101 (1-w 1 )e -α 1 τ 1 .(2)The system starts from (w 2 , y 0 ), goes through a τ 2 time forgetting condition and reaches (x 1 , y 1 ), then goes back to (w 2 , y 2 ) through a training condition.Here y 2 = M M +N ( 1 y , where M = (1+( e β 1 τ 2

Fig. 1 .

 1 Fig. 1. The ultimate stable states under the perturbation of measurement error


Christos P Loizou

Efthyvoulos C Kyriacou

Ioannis Seimenis

Marios Pantziaris

Christodoulos Christodoulou

Constantinos S Pattichis

Brain White Matter Lesions Classification in Multiple Sclerosis Subjects for the Prognosis of Future Disability

Keywords: MRI, multiple sclerosis, texture classification

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

Introduction

Multiple sclerosis (MS) is the most common autoimmune disease of the central nervous system, with complex pathophysiology, including inflammation, demyelination, axonal degeneration, and neuronal loss. Within individuals, the clinical manifestations are unpredictable, 1 particularly with regard to the development of disability [START_REF] Fazekas | The contribution of magnetic resonance imaging to the diagnosis of multiple sclerosis[END_REF].

Diagnostic evaluation of MS, performed by a specialized neurologist, is generally based on conventional magnetic resonance imaging (MRI) following the McDonald criteria [START_REF] Mcdonald | Recommended diagnostic criteria for multiple sclerosis: guidelines from the international panel on the diagnosis of multiple sclerosis[END_REF], and on clinical signs and symptoms. The development of modern imaging techniques for the early detection of brain inflammation and the characterization of tissue-specific injury is an important objective in MS research. Recent MRI studies have shown that brain and focal lesion volume measures, magnetization transfer ratio and diffusion weighted imaging-derived parameters can provide new information in diagnosing MS [START_REF] Bakshi | MRI in multiple sclerosis: current status and future prospects[END_REF]. Texture features quantify macroscopic lesions and also the microscopic abnormalities that may be undetectable using conventional measures of lesion volume and number [START_REF] Fazekas | The contribution of magnetic resonance imaging to the diagnosis of multiple sclerosis[END_REF], [START_REF] Kassner | Texture analysis: A review of neurologic MR imaging applications[END_REF]. Several studies have been published, where it was documented that texture features can be used for the assessment of MS lesions in: (i) differentiating between lesions for normal white matter (NWM), and the so called normal appearing white matter (NAWM) [START_REF] Harrison | MRI Texture analysis in multiple sclerosis: Toward a clinical analysis protocol[END_REF][START_REF] Herlidou-Meme | MRI texture analysis on texture test objects, normal brain and intracranial tumors[END_REF][START_REF] Mathias | Texture analysis of spinal cord pathology in multiple sclerosis[END_REF][START_REF] Loizou | Quantitative analysis of brain white matter lesions in multiple sclerosis subjects: Preliminary findings[END_REF][START_REF] Loizou | Quantitative analysis of brain white matter lesions in multiple sclerosis subjects[END_REF], and (ii) monitoring the progression of the disease over longitudinal scans [START_REF] Collewet | Influence of MRI acquisition protocols and image intensity normalization methods on texture classification[END_REF][START_REF] Yu | Distinct patterns of active and non-active plaques using texture analysis of brain NMR images in multiple sclerosis patients: Preliminary results[END_REF][START_REF] Zhang | Feature reduction and texture classification in MRI-Texture analysis of multiple sclerosis[END_REF][START_REF] Meier | Time-series analysis of MRI intensity patterns in multiple sclerosis[END_REF][START_REF] Loizou | Multi-scale Amplitude Modulation-Frequency Modulation (AM-FM) texture analysis of multiple sclerosis in brain MRI images[END_REF].

Our primary objective in this study is to develop texture classification methods that can be used to predict MS brain lesions that at a later stage are associated with advanced clinical disability and more specifically with the extended disability status scale (EDSS). Since the use of quantitative MRI analysis as a surrogate outcome is used as a surrogate measure in clinical trials, we hypothesise that there is a close relationship between the change in the extracted features and the clinical status and the rate of development of disability. We analyzed patient's images acquired at the initial stages of the disease-clinical isolated syndrome (CIS) (0 months) and we correlated texture findings with disability assessment scales. We interrelate therefore the EDSS scores [START_REF] Thompson | Multiple sclerosis: assessment of disability and disability scales[END_REF] with standard shape and texture features.

Preliminary findings of this study, for the texture analysis of NAWM in MS subjects were also published in [START_REF] Loizou | Quantitative analysis of brain white matter lesions in multiple sclerosis subjects: Preliminary findings[END_REF], [START_REF] Loizou | Quantitative analysis of brain white matter lesions in multiple sclerosis subjects[END_REF]. It should be note that the same problem was also investigated in a recent study published by our group using Amplitude Modulation-Frequency Modulation (AM-FM) analysis [START_REF] Loizou | Multi-scale Amplitude Modulation-Frequency Modulation (AM-FM) texture analysis of multiple sclerosis in brain MRI images[END_REF]. The motivation of this study was to investigate the usefulness of classical texture analysis, as well as compare the findings with the AM-M analysis.

MRI-based texture analysis was shown to be effective in classifying MS lesions from NWM and the so called, NAWM, with an accuracy of 96%-100% 5. In [START_REF] Herlidou-Meme | MRI texture analysis on texture test objects, normal brain and intracranial tumors[END_REF], the authors showed that texture features can reveal discriminant features for differentiating between normal and abnormal tissue, and for image segmentation. Significant differences in texture between normal and diseased spinal cord in MS subjects were found in [START_REF] Mathias | Texture analysis of spinal cord pathology in multiple sclerosis[END_REF] as well a significant correlation between texture features and disability. The median value increase of these texture features suggests that the lesions texture in MS subjects is less homogeneous and more complex than the corresponding healthy tissue (NWM) in healthy subjects [START_REF] Mathias | Texture analysis of spinal cord pathology in multiple sclerosis[END_REF]. Similar findings were also found in [START_REF] Loizou | Quantitative analysis of brain white matter lesions in multiple sclerosis subjects: Preliminary findings[END_REF] and [START_REF] Loizou | Quantitative analysis of brain white matter lesions in multiple sclerosis subjects[END_REF], where it was shown that median values of lesion texture features such as standard deviation, median, sum of squares variance, contrast, sum average, sum variance, difference variance, and difference entropy increased significantly with the progression of the MS disease when compared to NWM tissue. Statistical analysis (using the spatial gray level dependence matrices) has shown that there is a significant difference between lesions and NAWM or NWM. These findings may be beneficial in the research of early diagnosis and treatment monitoring in MS.

The differentiation between active and non-active brain lesions in MS subjects from brain MRI was also investigated in [START_REF] Yu | Distinct patterns of active and non-active plaques using texture analysis of brain NMR images in multiple sclerosis patients: Preliminary results[END_REF]. Here, it was shown that active lesions could be identified without frequent gadolinium injections, using run length analysis criteria. In [START_REF] Zhang | Feature reduction and texture classification in MRI-Texture analysis of multiple sclerosis[END_REF] the performance of texture analysis concerning discrimination between MS lesions, NAWM and NWM from healthy controls was investigated by using linear discriminant analysis. The results suggested that texture features can support early diagnosis in MS. When a combined set of texture features was used [START_REF] Zhang | Feature reduction and texture classification in MRI-Texture analysis of multiple sclerosis[END_REF], similar findings were reported.

In [START_REF] Meier | Time-series analysis of MRI intensity patterns in multiple sclerosis[END_REF], a pattern recognition system was developed for the discrimination of multiple sclerosis from cerebral microangiopathy lesions based on computer-assisted texture analysis of magnetic resonance images. It was shown that, MS regions were darker, of higher contrast, less homogeneous and rougher as compared to cerebral microangiopathy. Finally, in [START_REF] Loizou | Multi-scale Amplitude Modulation-Frequency Modulation (AM-FM) texture analysis of multiple sclerosis in brain MRI images[END_REF], the use of multi-scale Amplitude Modulation-Frequency Modulation (AM-FM) texture analysis of MS using magnetic resonance images from brain was introduced.

Materials & Methods

Material and MRI acquisition

Thirty subjects (15 males, and 15 females), aged 31.4±12.6 (mean age ± standard deviation) were scanned at 1.5T within one month following initial clinical evaluation to confirm CIS diagnosis. The transverse MR images used for analysis were obtained using a T2-weighted turbo spin echo pulse sequence (repetition time=4408 ms, echo time=100 ms, echo spacing=10.8 ms). The reconstructed image had a slice thickness of 5 mm and a field of view of 230 mm with a pixel resolution of 2.226 pixels per mm. Standardized planning procedures were followed during each MRI examination. The MRI protocol and the acquisition parameters were given in detail in [START_REF] Loizou | Quantitative analysis of brain white matter lesions in multiple sclerosis subjects: Preliminary findings[END_REF], [START_REF] Loizou | Quantitative analysis of brain white matter lesions in multiple sclerosis subjects[END_REF].

Initial clinical evaluation was made by an experienced MS neurologist (co author M. P.) who referred subjects for a baseline MRI upon diagnosis and clinically followed all subjects for over two years. All subjects were subjected to an EDSS test two years after initial diagnosis to quantify disability [START_REF] Thompson | Multiple sclerosis: assessment of disability and disability scales[END_REF].

A normalization algorithm was used to match the image brightness between the first (baseline) and the follow-up images (see [START_REF] Loizou | MRI intensity normalization in brain multiple sclerosis subjects[END_REF] for details). For this purpose, the neurologist manually segmented cerebrovascular fluid (CSF) areas as well as areas with air (sinuses) from all MS brain scans. Similarly, ROIs representing NWM, CSF and air from the sinuses were arbitrarily segmented from the brain scans of the 20 healthy subjects. The original image histogram was stretched, and shifted in order to cover all the gray scale levels in the image.

To introduce the objective of our study, an example in Fig. 1 is presented. Here, we show a transaxial T2-weighted MR image of a female patient at the age of 32, with an EDSS [START_REF] Thompson | Multiple sclerosis: assessment of disability and disability scales[END_REF] equal to 3. The image in Fig. 1a) corresponds to the initial diagnosis of a CIS of MS and the delineated lesion corresponding to the MS plaque is also depicted. Figure 1b) shows the magnified segmented lesion from a) (original images were acquired at a sampling rate of 2.226 pixels per mm). In what follows, texture analysis refers to the image processing of the extracted regions of interest (ROIs).

Manual delineations and visual perception

All MRI-detectable brain lesions were identified and segmented by an experienced MS neurologist and confirmed by a radiologist. Only well-defined areas of hyperintensity on T2-weighted MR images were considered as MS lesions. The neurologist manually delineated (using the mouse) the brain lesions by selecting consecutive points at the visually defined borders between the lesions and the adjacent NAWM on the acquired transverse T2-weigted sections. Similar regions corresponding to NAWM were delineated contralaterally to the detected MS lesions. The manual delineations were performed using a graphical user interface implemented in MATLAB developed by our group. Manual segmentation by the MS expert was performed in a blinded manner, without the possibility of identifying the subject or the clinical findings. The selected points and delineations were saved to be used for texture analysis.

Texture analysis

In this study texture features and shape parameters were extracted from all MS lesions detected and segmented [START_REF] Loizou | Quantitative analysis of brain white matter lesions in multiple sclerosis subjects: Preliminary findings[END_REF], [START_REF] Loizou | Quantitative analysis of brain white matter lesions in multiple sclerosis subjects[END_REF], while texture features were also calculated for the NAWM ROIs. The overall texture and shape features for each subject were then estimated by averaging the corresponding values for all lesions for each subject. The following texture features were extracted: (i) Statistical Features: a) mean, b) variance, c) median value, d) skewness, e) kurtosis, f) energy and g) entropy. (ii) Spatial Gray Level Dependence Matrices (SGLDM) as proposed by Haralick et al.

[17]: a) angular second moment (ASM), b) contrast, c) correlation, d) sum of squares variance (SOSV), e) inverse difference moment (IDM), f) sum average (SA), g) sum variance (SV), h) sum entropy (SE), i) entropy, j) difference variance (DV), k) difference entropy (DE), and l) information measures of correlation (IMC). For a chosen distance d (in this work d=1 was used) and for angles θ = 0 0 , 45 0 , 90 0 , and 135 0 , we computed four values for each of the above texture measures. Each feature was computed using a distance of one pixel. Then for each feature the mean values and the range of values were computed, and were used as two different features sets.

(iii) Gray Level Difference Statistics (GLDS) [START_REF] Weszka | A comparative study of texture measures for terrain classification[END_REF]: a) homogeneity, b) contrast, c) energy, d) entropy, and e) mean. The above features were calculated for displacements δ=(0, 1), ( (vi) Laws Texture Energy Measures (LTEM) [START_REF] Weszka | A comparative study of texture measures for terrain classification[END_REF]: LL-texture energy from LL kernel, EE-texture energy from EE-kernel, SS-texture energy from SS-kernel, LE-average texture energy from LE and EL kernels, ES-average texture energy from ES and SE kernels, and LS-average texture energy from LS and SL kernels. (vii) Fractal Dimension Texture Analysis (FDTA) [START_REF] Mandelbrot | The Fractal Geometry of Nature[END_REF]: The Hurst coefficients for dimensions 4, 3 and 2 were computed. (viii) Fourier Power Spectrum (FPS) [START_REF] Weszka | A comparative study of texture measures for terrain classification[END_REF], [START_REF] Wu | Texture features for classification of ultrasonic images[END_REF]: a) radial sum, and b) angular sum. (vii) Shape Parameters: a) X-coordinate maximum length, b) Y-coordinate maximum length, c) area, d) perimeter, e) perimeter2/area, f) eccentricity, g) equivalence diameter, h) major axis length, i) minor axis length, j) centroid, k) convex area, and l) orientation.

Models support vector machines classification

Classification analysis was carried out to classify brain MS lesions delineated on the baseline MRI scans into two classes according to the EDSS score that each patient was allocated two years following initial diagnosis: (i) MS subjects with EDSS≤2, and (ii) MS subjects with an EDSS>2. Thus, the classification goal was to differentiate between lesions that were subsequently associated with mild (EDSS≤2) or advanced disability (EDSS>2). The classification was applied on 30 subjects (15 males and 15 females).

Each classifier was implemented in Matlab, using Support Vector Machines (SVM). The SVM network was investigated using Gaussian Radial Basis Function kernels [START_REF] Cristianini | An Introduction to Support Vector Machines and Other Kernel-based Learning Methods[END_REF]; this was decided as the rest of the kernel functions could not achieve satisfactory results. The Parameters for the classifier were selected using ten-fold cross validation for each set. Parameter C of SVM was 1. The leave-one-out estimate was used for validating the classification models.

The performance of the classifier models were measured using the receiver operating characteristics (ROC) curve parameters [START_REF] Ebrchart | Neural Networks PC Tools A Practical Guide[END_REF]. The parameters calculated for ROC curves are: (i) true positives (TP) when the system correctly classifies subjects as EDSS≤2, and EDSS >2, (ii) false positives (FP) where the system wrongly classifies subjects as EDSS>2 while they are in the group of EDSS ≤2, (iii) true negatives (TN) when the system correctly classifies subjects as EDSS ≤2. We also compute the Sensitivity (SE) which is the likelihood that a subject with EDSS >2 will be detected given that it is EDSS>2 and Specificity (SP) which is the likelihood that a subject will be classified as EDSS ≤2 given that he/she is EDSS ≤2. For the overall performance, we provide the correct classification (CC) rate which gives the percentage of correctly classified subjects.

Results

We present in Fig. 2 box plots of the lesion texture features mean (Mean0_1, Mean0_2) median (Median0_1, Median0_2) and sum average (SA0_1, SA0_2) for 0 months for EDSS≤2 (_1) (N=18) and EDSS>2 (_2) (N=11). Statistical analysis was carried out for all texture features given in section 2.3 for lesions and NAWM for EDSS≤2 versus EDSS>2. The results of selected features, including those that showed significant difference are tabulated in Table 1. For lesions, features mean, median, IDM and SA could differentiate between EDSS≤2 and EDSS>2 (see also Fig. 2). For the NAWM it is shown that there are no significant differences between subjects at 0 months with EDSS>2 and EDSS≤2, for all texture features.

Table 2 presents the lesion classification results. According to Table II the best correct classification rate was achieved as a combination of GLDS, NGTDM and SFM feature sets and was 77%.

Discussion

The objective of our study was to investigate texture classification analysis in an effort to differentiate between MS brain lesions that, in a subsequent stage, will be associated with advanced diseased disability from those lesions that will be related to a mild disability score.

Table 1 showed that there are some texture features that can be possibly used to differentiate between subjects with mild (EDSS≤2) and advanced (EDSS>2) MS disease states. These features are mean, median, IDM and SA. Table 2 showed that the combination of GLDS, NGTDM & SFM texture features gave a % of CC rate of 77% when classifying subjects according to their long-term disability status.

These results can be compared with another study carried out by our group [START_REF] Loizou | Multi-scale Amplitude Modulation-Frequency Modulation (AM-FM) texture analysis of multiple sclerosis in brain MRI images[END_REF], where AM-FM feature analysis was used but on a slightly larger set of subjects. The findings suggest that AM-FM characteristics succeed in differentiating (i) between NWM and lesions, (ii) between NAWM and lesions, and (iii) between NWM and NAWM. A support vector machine (SVM) classifier succeeded in differentiating between patients that, two years after the initial MRI scan, acquired an EDSS ≤ 2 from those with EDSS > 2 (correct classification rate = 86%). The best classification results were obtained from including the combination of the low-scale instantaneous amplitude (IA) and instantaneous frequency (IF) magnitude with the medium-scale IA. Among all AM-FM features, the medium frequency IA histogram alone gave the best AUC results (AUC=0.76). Also the best % of CC achieved was 86%. The IA histograms from all frequency scales contributed to a multi-classifier system that gave 86% correct classification rate. The best classifier results also used the IF magnitude from both the low and the high frequency scales.

To the best of our knowledge, no other studies were carried out for differentiating between the aforementioned two disability scores. Several studies were carried out for differentiating and classifying NWM, and or NAWM, and lesions, as these are summarized below.

Various studies have been performed in order to establish a relationship between the various gray levels and texture features [START_REF] Kassner | Texture analysis: A review of neurologic MR imaging applications[END_REF]- [START_REF] Yu | Distinct patterns of active and non-active plaques using texture analysis of brain NMR images in multiple sclerosis patients: Preliminary results[END_REF]. In [START_REF] Harrison | MRI Texture analysis in multiple sclerosis: Toward a clinical analysis protocol[END_REF] MRI texture analysis based on statistical, autoregressive model, and wavelet-derived texture analysis was performed on 30 MS subjects. The classification accuracy between MS lesions, NAWM and plaques NWM, was 96%-100%.

Likewise in [START_REF] Zhang | Texture analysis of multiple sclerosis: a comparative study[END_REF], texture analysis was performed on MR images of MS subjects and normal controls and a combined set of texture features were explored in order to better discriminate tissues between MS lesions, NAWM and NWM. The results demonstrated that with the combined set of texture features, classification was perfect (100%) between MS lesions and NAWM (or NWM), less successful (88.89%) among the three tissue types and worst (58.33%) between NAWM and NWM. Furthermore, it was shown that compared with GLCM-based features, the combined set of texture features were better at discriminating MS lesions and NWM, equally good at discriminating MS lesions and NAWM and at all three tissue types, but less effective in classification between NAWM and NWM. The study suggested that texture analysis with the combined set of texture features may be equally good or more advantageous than the commonly used GLCM-based features alone in discriminating MS lesions and NWM/NAWM and in supporting early diagnosis of MS.

Yu et al. [START_REF] Yu | Distinct patterns of active and non-active plaques using texture analysis of brain NMR images in multiple sclerosis patients: Preliminary results[END_REF] performed textural feature analysis to discriminate active and nonactive MS lesions in a study of 8 subjects with relapsing remitting MS, by using linear discriminant analysis using 42 first-and second order statistical textural features. Applications of the run-length method have been very limited compared with other methods, yet Yu et al. found that run-length matrix features actually outperformed gray-level co-occurrence matrix features in the identification of active MS lesions, with run-length matrix features distinguishing active from inactive lesions with 88% sensitivity and 96% specificity. Conversely, none of the gray-level co-occurrence features provided any discrimination between lesion subtypes.

In another study [START_REF] Zhang | Feature reduction and texture classification in MRI-Texture analysis of multiple sclerosis[END_REF], it was shown that texture analysis can achieve high classification accuracy (>=90%) in tissue discrimination between MS lesions and NAWM.

In [START_REF] Mayerhoefer | Are signal intensity and homogeneity useful parameters for distinguishing between benign and malignant soft tissue masses on MR images? Objective evaluation by means of texture analysis[END_REF] the authors tried to objectively identify possible differences in the signal characteristics of benign and malignant soft tissue masses on MRI by means of texture analysis and to determine the value of these differences for computer-assisted lesion classification. Fifty-eight subjects with histologically proven soft tissue masses (benign, n=30; malignant, n=28) were included and their masses were texture analyzed. The best results of soft tissue masses classification were achieved using texture information from short tau inversion recovery images, with an accuracy of 75.0% (sensitivity, 71.4%; specificity, 78.3%) for the k-NN classifier, and an accuracy of 90.5% (sensitivity, 91.1%; specificity, 90.0%) for the artificial neural network classifier. The authors concluded that, texture analysis revealed only small differences in the signal characteristics of benign and malignant soft tissue masses on routine MRI.

Future work

Further research work on a larger number of subjects is required for validating the results of this study and for finding additional shape and texture features that may provide information for longitudinal monitoring of the lesions on the initial MRI scan of the brain of patients with CIS; and improve the final classification rate. Additionally the use of other classifiers could suggest a better classification a scheme and will be examined. Finally the proposed methodology could be possibly used for the assessment of subjects at risk of developing future neurological events and disease progression as measured by increased EDSS score.

Fig. 1 .

 1 Fig. 1. a) ROI drawn on MR image of the brain obtained from a 32 year old female MS patient with an EDSS=3 at 0 months, and b) magnified segmented lesion that was acquired at a pixel resolution of 2.226 pixels per mm. The bar below the lesion shows the size of 10 mm. The grey scale median and inter-quartile range (IQR) of the segmented lesion were 108 and 9.6, respectively.

Fig. 2 .

 2 Fig. 2. Box plots for the median ± interquartile range (IQR) values for texture features mean, median and sum average, from MS lesions at 0 months for EDSS<=2 (N=18) and EDSS>2 (N=11) corresponding to feature notation _1, and _2 respectively (see also TABLE I). Interquartile range (IQR) values are shown above the box plots. In each plot we display the median, lower, and upper quartiles and confidence interval around the median. Straight lines connect the nearest observations within 1.5 of the IQR of the lower and upper quartiles.

     , and their mean values were taken. (iv) Neighbourhood Gray Tone Difference Matrix (NGTDM) [19]: a) coarseness, b) contrast, c) busyness, d) complexity, and e) strength. (v) Statistical Feature Matrix (SFM) 20: a) coarseness, b) contrast, c) periodicity, and d) roughness.

Table 1 .

 1 Statistical Analysis Of The Texture Features For Lesions And NAWM Subjects With EDSS≤2 (N=18) and EDSS>2 (N=11) Based on The Mann Whitney (first row) and the Wilcoxon (2nd row) Rank Sum Tests at p<0.05. Significant Difference is Depicted With S And Non Significant Difference is Depicted With NS. The p Values Are Shown in Parentheses.

	Feature	Lesions	NAWM
	Mean	S (0.03)	NS (0.18)
	Median	S (0.03)	NS (0.17)
	STD	NS (0.83)	NS (0.07)
	Contr	NS (0.16)	NS (0.36)
	SOSV	NS (0.83)	NS (0.55)
	Entr	NS (0.06)	NS (0.58)
	IDM	S (0.01)	NS (0.36)
	SA	S (0.04)	NS (0.17)
	SV	NS (0.83)	NS (0.19)
	DV	NS (0.44)	NS (0.18)
	DE	NS (0.35)	NS (0.77)
	NAWM: Normal appearing white matter, STD: Standard deviation, Contr.: Contrast, SOSV: Sum of
	square variance, Entr: Entropy, IDM: Inverse difference moment, SA: Sum average, SV: Sum
	variance, DE: Difference entropy.

Table 2 .

 2 Classification Results Using The Texture And Shape Feature Set And The Support Vector Machines Classifier For Two Classes: EDSS ≤2 and EDSS>2 (N=30). The Results Are Presented Using the ROC Measures: Percentage of correct classifications (%CC), percentage of false positives (%FP), percentage of false negatives (%FN), percentage sensitivity (%SE) and percentage specificity (%SP).

	Feature Group	%CC	%FP	%FN	%SE	%SP
	FOS	43	47	67	33	53
	SGLDM(mean)	63	53	20	80	47
	SGLDM(range)	63	53	20	80	47
	GLDS	67	47	20	80	53
	NGTDM	67	40	27	73	60
	SFM	67	40	27	73	60
	TEM	67	40	27	73	60
	FDTA	60	33	7	53	67
	FPS	60	53	27	73	47
	NGTDM & SFM	73	27	27	73	73
	GLDS & NGTDM & SFM	77	27	20	80	73
	SHAPE	53	53	40	60	47

6 Acknowledgment This work was supported through the project Quantitative and Qualitative Analysis of MRI Brain Images ΤΠΕ/ΟΡΙΖΟ/0308(ΒΙΕ)/15, 12/2008-12/2010, of the Program for Research and Technological Development 2007-2013, co funded by the Research Promotion Foundation of Cyprus.


Jaime R Ticay-Rivas

Marcos Del Pozo-Baños

Carlos M Travieso
email: ctravieso@idetic.eu

Jorge Arroyo-Hernández

Santiago T Pérez
email: sperez@idetic.eu

Jesús B Alonso
email: jalonso@idetic.eu

Federico Mora-Mora

Marcos Del Pozo-Baños

Jorge Arroyo- Hernández

Pollen Classification Based on Geometrical, Descriptors and Colour Features Using Decorrelation Stretching Method

Keywords: Pollen grains, pollen classification, colour features, geometrical features, neural networks, decorrelation stretch

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

Introduction

Over 20% of all the world's plants are already at the edge of becoming extinct [START_REF]Plants under pressure: a global assessment[END_REF]. Saving earth's biodiversity for future generations is an important global task [START_REF] Sytnik | Preservation of biological diversity: Top-priority tasks of society and state[END_REF] and many methods must be combined to achieve this goal. Saving flora biodiversity involves mapping plant distribution by collecting pollen and identifying them in a laboratory environment. Pollen grain classification is a qualitative process, involving observation and discrimination of features.

The manual method, heavily depending on experts, is still the most accurate and effective, but takes considerable amounts of time and resources, limiting research progress [START_REF] Stillman | The Needs and Prospects for Automation in Palynology[END_REF]. Therefore, automatic recognition of pollen species by means of computer vision is a highly prioritized subject in palynology.

A lot of progress has been made in the field of automatic pollen species classification. In [START_REF] Li | Pollen texture identification using neural networks[END_REF], it showed one of the first works in which texture features and neural networks were used for pollen grains identification task. In [START_REF] France | A new approach to automated pollen analysis[END_REF] detection and classification of pollen grains were identified as the main and difficult tasks. In both cases neural networks were used. The characteristics of pollen grains for their correct classification were studied in different studies. In [START_REF] Treloar | Towards Automation of Palynology 1: Analysis of Pollen Shape and Ornamentation using Simple Geometric Measures, Derived from Scanning Electron Microscope Images[END_REF], the shape and ornamentation of the grains were analyzed, using simple geometric measures. Concurrence matrices were applied for the measurement of texture in [START_REF] Li | Towards Automation of Palynology 2: The Use of Texture Measures and Neural Network Analysis for Automated Identification of Optical Images of Pollen Grains[END_REF]. Both techniques were combined in [START_REF] Zhang | Towards Automation of Palynology 3: Pollen Pattern Recognition using Gabor Transforms and Digital Moments[END_REF] and [START_REF] Rodriguez-Damian | Automatic detection and classification of grains of pollen based on shape and texture[END_REF]. In [START_REF] Rodriguez-Damian | Pollen classification using brightness-based and shape-based descriptors[END_REF] brightness and shape descriptors were used as vector features.

Another interesting approach is found in [START_REF] Boucher | Object recognition from 3D blurred images[END_REF]. Here, a combination of statistical reasoning, feature learning and expertise knowledge from the application was used. Feature extraction was applied alternating 2D and 3D representations. Iterative refinement of hypotheses was used during the classification process. A full 3D volume database of pollen grain was recorded in [START_REF] Ronneberger | General-purpose object recognition in 3D volume data sets using gray-scale invariants -classification of airborne pollen-grains recorded with a confocal laser scanning microscope[END_REF] using a confocal laser scanning microscope. 14 invariant gray-scale features based on integration over the 3D Euclidian transformation were extracted and used for classification. A more recent work with an ambitious approach is shown in [START_REF] Allen | Machine vision for automated optical recognition and classification of pollen grains or other singulated microscopic objects[END_REF]. It describes an automatic optical recognition and classification of pollen grains system. This is able to locate pollen grains on slides, focus and photograph them, and finally identify the species applying a trained neural network.

There is no doubt that significant progress has been made, but more knowledge is needed in the combination of features that have been studied and colour features in order to increase the number of species that systems are capable to recognize and classify with a high rate of success. This work introduces the use of colour features over decorrelated images. A novel combination of these features with geometrical measures, and Fourier descriptors of morphological details are proposed as pollen grains discriminative features, being this our first work in this area. The remainder of this paper is organized as follows. In Section 2 the process of pollen grains detection is explained. Section 3 introduces the techniques used for feature extraction. The classification block is presented in Section 4. Experimental methodology and its results are exposed in Section 5. Finally in Section 6, discussions and conclusions are shown.

Pollen Grains Detection

Since there is not a uniform focus and conglomerates may appear in the images, a simple semi-automatic algorithm was developed to capture the pollen grains, indicating and cropping the pollen image. Thus, only one grain appears in each pollen grain sub-image. An example of each pollen class and its Scientist name (family, gender and specie by column) can be seen in figure 1.

Once sub-images of the pollen grains have been selected, automatic detection is performed. First, the colour contrast of the image is increased in order to locate the pollen grain in the image. Then a decorrelation stretching (DS) of the image is done. This technique is explained in detail in [START_REF] Cambell | The decorrelation stretch transformation[END_REF] and [START_REF] Amanatiadis | A comparative study of invariant descriptors for shape retrieval[END_REF]. The main idea is to achieve a variance-equalizing enhancement in order to uncorrelate and to fill the colour space.

The pollen grain image is represented by a matrix with dimensions M×N×3 where M×N represents the image dimensions of each colour channel. Thus, the pixels of the image are represented by the vector

I x   . Let V express
the covariance or correlation matrix of x and the eigenanalysis associated [START_REF] Amanatiadis | A comparative study of invariant descriptors for shape retrieval[END_REF] with the principal component analysis described as:

t UEU V  (1)
Where the columns of U denote the eigenvectors and the diagonal elements of the diagonal matrix E denote the corresponding eigenvalues. Using the transform t U the vector x is transformed in a new vector space where the new data y is given by:

x U y t  (2)
The scaling, applied in order to equalize the variance, is achieved in this space dividing each value of y by its corresponding standard deviation The transformation to the original space is given by

x U UE z t 2 / 1   (5)
The resulting image is transformed to the HSV colour model (Hue, Saturation and Value). The S channel is used in order to take advantage of the difference of brightness between the background and the pollen grain. Histogram equalized image is necessary before the binarization of the image. The next factor is multiplied to the histogram equalization:

8 . 0 / 1 ) 5 . 1 )* 1 . 0 ((    m e f ( 6
)
where m is the mean of the S channel insuring that the equalization is made in a correct range.

The binarization is applied setting the threshold at 0.45. Morphological methods like close, erosion, dilation, etc. are used in order to consider only the area where the pollen grain is found. Figure 2 shows the obtained images:

Discriminative Feature Extractions

Discriminative feature extractions were divided in three steps and are described in the next sections:

Geometrical Features

The geometrical features have been used in almost all works related to grain classification. The area, convex area, and perimeter are used in this work as they have been proven to provide a fair amount of the basic information of the pollen grain. Area is the number of pixels in the considered region. Convex area is the number of pixels of the binary image that specifies the convex hull. Perimeter is calculated as the distance between each adjoining pair of pixels around the border of the pollen grain.

Fourier Descriptors

The Fourier Descriptor represents the shape of the pollen grain in the frequency domain [START_REF] Rodriguez-Damian | Pollen classification using brightness-based and shape-based descriptors[END_REF], [START_REF] Amanatiadis | A comparative study of invariant descriptors for shape retrieval[END_REF]. Each boundary point is represented as a complex number. In order to eliminate the bias effect, the centroid of the shape of the pollen grain is introduced in the complex representation:

] [ ] [ c c y y i x x s     (7)
The contour is sampled each 2 degrees using the Cartesian to radian coordinates transformation. The Fourier Descriptor is defined as the power spectrum of the discrete Fourier Transform of s. In this work the Discrete Cosine Transformation (DCT) is applied to the resulting power spectrum in order to select the most significant values of the Fourier Descriptor. Finally, the first three coefficients of the transformation were used as features on this approach.

Colour Features

In section 2, decorrelation stretching method (DS) is applied to the pollen grain images. The averages of R (Red), B (Blue), H (Hue) channels were calculated as three colour parameters. It is important to mention that only the pollen grain region was considered. The Gray Level Co-occurrence Matrix (GLCM) from H channel was used in order to obtain its contrast. The distance and direction setting in the GLCM calculation was 1 pixel and 0 degree respectively.

Classification System

Artificial neural networks (NN) have been used in several studies [START_REF] Li | Pollen texture identification using neural networks[END_REF], [START_REF] France | A new approach to automated pollen analysis[END_REF], [START_REF] Li | Towards Automation of Palynology 2: The Use of Texture Measures and Neural Network Analysis for Automated Identification of Optical Images of Pollen Grains[END_REF], [START_REF] Zhang | Towards Automation of Palynology 3: Pollen Pattern Recognition using Gabor Transforms and Digital Moments[END_REF], [START_REF] Allen | Machine vision for automated optical recognition and classification of pollen grains or other singulated microscopic objects[END_REF] demonstrating good results. Multi-Layer Perceptron NN trained by Back Propagation algorithm (MLP-BP) were used. Neurons implemented the hyperbolic tangent sigmoid transfer function. The number of neurons in the input layer is the same as the number of selected features. Following experimentation, it was determined to use one hidden layer. The number of neurons that showed better results ranged between 20 and 80. The best results were found using 50 neurons. The numbers of neurons in the output layer was 17, that is, the number of classes. 5000 iterations were performed in the training phase.

Since weights of neurons in the input layer were initialized randomly, 30 NNs were trained in order to determine the classes. Two decision techniques (DT) were applied: class majority voting (CMV) and adding scores (AS). The former selects the most repeated result along all 30 NNs. These results at every NN are computed as the index of the maximum scoring neuron. The AS adds the output scores of corresponding neurons from each 30 NNs, before selecting the index of the maximum value.

Experiments and Results

Database

The database contains a total of 426 pollen grain sub-images as those presented in figure 1. These sub-images correspond to 17 genders and species from 11 different families of tropical honey plants from Costa Rica, Central America.

Experimental Methodology and Results

Parameters introduced in section 3 were tested for both raw images and DS processed images. For each of them, apart from the original dataset of parameters (ORIG), a second dataset of parameters (ORIGPCA) was obtained by applying PCA to the original dataset. About the classification, CMV and AS decision techniques were applied separately. Moreover, both parameter sets were fused at the score level using again the AS technique. This makes a total of 8 combinations of raw/DS-processed images, original/PCA data set, and CVM/AS decision technique.

In order to obtain valid results, a series of cross-validation methods were performed on every experiment: 30 iterations of 50% hold out (50-HO), and k-folds with k equal 3, 4, and 5 (3-, 4-, and 5-folds). A second set of features were obtained applying the DS method to the original images. Tables 1 and2 show all obtained results for both original and DS processed set of images, showing results as mean (%) and standard desviation (±) from 30 iterations.

Param. DS DT 50-HO 3-folds 4-folds 5-folds ORIG + ORIGPCA

No AS 91.19%  2.07 92.96%  1.74 94.61%  3.05 95.78%  4.45 ORIG + ORIGPCA Yes AS 95.33%  1.31 96.53%  0.21 96.63%  1.17 97.48%  1.97

Discussion and Conclusions

This work has introduced the use of colour features over decorrelated stretched images, for enhance pollen grain classification. More standard characteristics like geometrical features and Fourier descriptors have been added to the pollen grain descriptions. Over this multiple feature vector, PCA has proven to increase the system's performance.

Experimental results shown in Table 1 indicate how the use of DS technique significantly improves the recognition of pollen grain species. Moreover, this improvement can be seen in terms of success rate and stability (standard deviation). This can be explained by the fact that colours channel variances are equalized and uncorrelated, reducing intra-class variation and differentiating inter-class samples. Furthermore, the use of PCA on the original database increased the success rate. This is consistent with the fact that PCA transformation projects samples over more meaningful dimensions, maximizing variance.

Finally, both data sets (ORIG and ORIGPCA) were fused using the AS score fusion, obtaining a maximum success rate of 96.4%, and a standard deviation of 1.16 with 17 pollen species. It is worth to mention that these results improve those achieved by other authors [START_REF] Li | Towards Automation of Palynology 2: The Use of Texture Measures and Neural Network Analysis for Automated Identification of Optical Images of Pollen Grains[END_REF] [8] [9] [START_REF] Rodriguez-Damian | Pollen classification using brightness-based and shape-based descriptors[END_REF], even though the number of classified species was significantly larger.

Fig. 1 .

 1 Fig. 1. Examples of each pollen class.

Fig. 2 .

 2 Fig. 2. Images obtained in the pollen grains detection.

Table 1 .

 1 Results obtained for the original set of parameters and after applying PCA, for both CMV and AS score fusion methods.No No CMV 86.25%  2.67 87.98%  3.13 90.01%  3.77 90.44% 3.53 No Yes CMV 86.13%  3.00 88.98%  7.54 91.79%  7.01 92.62%  7.17 No No AS 86.48%  2.75 88.87%  2.70 90.30%  3.45 90.19%  3.84 No Yes AS 86.55%  3.08 90.25%  3.36 93.34%  4.74 93.80%  4.25 Yes No CMV 88.08%  1.73 90.91%  1.11 90.98%  2.45 91.37% 3.34 Yes Yes CMV 93.99%  1.63 95.22%  0.24 94.22%  2.76 95.73% 1.84 Yes No AS 88.20 %  2.19 90.91%  0.62 90.86%  2.12 91.76%  3.27 Yes Yes AS 94.18 %  1.71 95.55%  0.15 95.36%  1.88 95.73% 1.84

	DS PCA DT 50-HO	3-folds	4-folds	5-folds

Table 2 .

 2 Results obtained merging ORIG and ORIGPCA parameters using AS score fusion.

Acknowledgements. This work has been supported by Spanish Government, in particular by "Agencia Española de Cooperación Internacional para el Desarrollo" under funds from D/027406/09 for 2010, and D/033858/10 for 2011. To M.Sc. Luis A. Sánchez Chaves, Tropical Bee Research Center (CINAT) at Universidad Nacional de Costa Rica, for provides the database images and palynology knowledge.


Amin Hosseinian-Far

Elias Pimenidis

Hamid Jahankhani

D C Wijeyesekera

Financial Assessment of London Plan Policy 4A.2 by Probabilistic Inference and Influence Diagrams

Keywords: London Plan, Climate Change Mitigation, Influence Diagrams, Probabilistic Inference, Analytica, Sustainability

London Plan is the London mayor's Spatial Development Strategy. This strategic long-term plan comprises of proposals for different aspects of change within the London boundary. Furthermore, the proposals include chapters outlining adjustments in each facet. Policy 4A.2 reflects the Climate Change Mitigation scheme. Some consultations and research works have been performed to this point, but an extensive cost assessment has not been done. This paper reflects a financial assessment by means of Influence Diagrams based upon the London Plan policies 4A.X.

Introduction

There are various definitions for Sustainable City which is also known as eco-city. They generally denote to a city which tries to have the minimum energy, and food input and has the minimum CO2 emission, energy waste, general waste, methane and water pollution. The term was first introduced by Richard Register in 1987 as "An Ecocity is an ecologically healthy city" [START_REF] Register | Ecocity Berkeley, Building Cities for a Healthy Future[END_REF]. There are different plans for attaining the status of a sustainable city; they even vary from town to town and city to city. The United Kingdom is among the contributors to ecocity plans. The aim of this paper is to demonstrate how influence diagrams can be used to do knowledge representation for a strategic, uncertain scenario. The case study selected is the London Plan Policy 4A.2. Section 1 discusses the key features of the London Plan relevant to this work; Section 2 outlines Influence Diagrams and the advantages of using them for this scenario. Finally the introduction to ID implementation is considered.

London Plan

The London Plan is the Mayor's Spatial Development Strategy, which covers different aspects of London's development; and a draft version of the plan has been out for public consultation since 2006 [START_REF]The Mayor's London Plan Briefing Report. Friends of Earth[END_REF]. The section of the plan that this research focuses on is that of tackling the climate change and addressing sustainability analysis for London as a city. It is a strategic plan using a planning system and a set of aims and objectives relating to Energy Strategy and the control of climate change. The policies in the London Plan outline the proposals and schemes from London Transport strategic arrangement to the city Design. Climate Change is a sub division of the overall proposal.

Climate Change Policy Area

Climate Change as a sub division of the London Plan includes policies concerning different aspects of climate change from Air, Noise, Tackling Climate Change and even to Hazardous Substances. The Tackling Climate Change is the main concern of this research study. The policies involved are 4A.1 Tacking Climate Change, 4A.2 Mitigating Climate Change, and finally 4A.3 Sustainable Design and Construction. In spite of the categorization, these policies have overlapping areas between each other. The basis for the model which will be described in this work is Policy 4A.2, but some statistics and proposals will be borrowed from 4A.1 and 4A.3.

Tackling Climate Change

The main focus of this policy is on 60 % carbon dioxide reduction in London by 2050 relatively to the 1990 base. Comparatively, 15%, 20%, 25%, and 30% reduction targets have been set for 2010, 2015, 2020, and 2025 [START_REF]The Mayor's London Plan Briefing Report. Friends of Earth[END_REF]. According to the plan draft new Building Development regulations are needed to be introduced; not only the new built properties should be regulated, but also the existing houses and properties should adapt to the criteria accordingly. According to Day, Ogumka, & Jones 20% out of the 60% reduction can be gained through the use of renewable energies [START_REF] Day | A:The Use of the Planning System to Encourage Low Carbon Energy Technologiesin Buildings in[END_REF]. Foreman Roberts research debates that the 20% out of 60% solely by renewables is not achievable [START_REF] Roberts | Renewables and the London Plan[END_REF]. Despite of all discussions, the cost assessment will be performed by only considering the policies published by the Greater London Authority (GLA). There are 147 case studies by further plan analysis, as approved by the London mayor. One of the main concerns of the plan is that of the use of Combined Heat and Power (CHP). CHP utilization would the most cost effective approach in the London Plan [START_REF] Cibse | Small-scale combined heat and power for building[END_REF].

The other two policies focus on more thorough details of the London Climate Change route. According to research conducted by London South Bank University the impact assessment of the approach shown in Figure 1 should prove suitable to assess the general perspective of the London Plan climate change policies: Fig. 1. Climate Change policy in three main phase [START_REF] Day | Renewable Energy in the City[END_REF] Use less energy (be lean) supply energy efficiently (be clean)

Use renewable energy (be green)

According to the same research, the level of CO2 emission savings based on Energy Efficiency is higher than what could be expected by the use of renewable energy, as shown in table 1 below. One of the debates relating to the London Plan (LP) is that of the targets set by the government. The government has set 60% CO2 reduction by 2050; below with the specific figure of 60% being quite unclear as to what the government envisages. Fig. 2 illustrates the most predominant understanding and it this viewpoint that the work presented here adopts: Fig. 2. Calculation of energy/carbon dioxide savings [START_REF]The London Plan, Spatial Development Strategy for Greater London, Consolidated with Alterations since 2004[END_REF] Energy Efficiency is the first phase according to the plan and also based on the LSBU research data is shown as having the highest reduction. The use of Combined Heat and Power (CHP) which is also called cogeneration refers to generation of electricity as well as heating the development. Renewable energy systems do not play a major role in the reduction process, but the introduction of more efficient renewables could yield an advantage. The application of the plan on 260 approved developments has revealed the following results: Fig. 3. Savings and approaches [START_REF] Day | Renewable Energy in the City[END_REF] 3 Influence Diagrams An Influence Diagram (ID) is a network visualization for probabilistic models [START_REF] Shachter | Probabilistic Inference and Influence Diagrams[END_REF]. It also includes decision analysis and can be easily converted to decision trees [START_REF] Ronald | Influence Diagrams in: Decision Analysis[END_REF]. Howard and Matheson initially introduced the concept in 1981 and according to them Influence Diagrams can be considered as a bridge between qualitative and quantitative analysis [START_REF] Howard | Influence diagrams[END_REF]. The concept of ID has had substantial influence on Artificial Intelligence (AI) and has made AI more mature [START_REF] Boutilier | The Influence of Influence Diagrams on Artificial Intelligence in: Decision Analysis[END_REF].

An influence diagram consists of nodes and arcs. The nodes and arcs form a visual network [START_REF] Shachter | Model Building with Belief Networks and Influence Diagram in: Advances in Decision Analysis: From Foundations to Applications by Edwards[END_REF]. The nodes may be of different types: Variables, Constants, Decision Nodes, Chance Node, and Objective; although there are other naming conventions for the ID nodes. They can also be classified as deterministic or probabilistic nodes. The Arcs specify the relationship and dependencies of the variables [START_REF]Influence Diagrams[END_REF]. Definition of each node determines an associated mathematical expression for each node.

Influence Diagrams are a proper tool for decision making and analyzing a scenario under uncertainty [START_REF] Liao | Efficient Non-myopic Value-of-information Computation for Influence Diagrams[END_REF]. Applying ID to London Plan policy analysis is an appropriate approach as the London Plan involves some uncertainty points. Furthermore Influence Diagrams can be furthermore broken down and illustrate more details. In complex situations such as the London Plan, the boundaries and level of details are in the hands of the analyst and can even be further refined later.

A Well-Formed Influence Diagram (WFID) is an ID which is not only well visualized, but also includes well-defined mathematical expressions for the nodes [START_REF] Sage | Systems Engineering[END_REF]. In addition WFIDs can be evaluated using reversal and removal approaches. This can be used for complex inferential and probabilistic scenarios [START_REF] Shachter | Evaluating Influence Diagrams in: Operations Research[END_REF].

In terms of categorization and classification of the subject, Influence Diagrams are classified under 'Bayesian Updating, Certainty Theory and Fuzzy Logic' category (Fig. 4).

Implementation

Indicators

A set of indicators for the London Plan monitoring phase is introduced, but the list does not cover the 4A.X policies. Therefore the indicators and variables used in the model are derived from the LP qualitative and quantitative statements within the boundary of this study.

Model Creation

The Integrated Development Environment (IDE) used for this research is called Analytica, many academic and industrial users of the software review Analytica as one of the most powerful applications for modeling using Influence Diagrams; even though other modeling environments can be used.

There are various recommendations for the starting point of influence diagrams. A systemic approach is to create the systems map first, and based on the systems map the influence diagram can be created. By using the systems map, the areas of the analysis and boundaries of the system under investigation will be determined [START_REF]Influence Diagrams[END_REF]. In addition the influential parts of the research would be easier to find. A tentative and general systems map for this research is outlined in Fig. 6: The systems map illustrated above draws a boundary for the system. The Environmental and Social assessments of the LD policy 4A.X overlap with these two domains. The economic or financial assessment domain is fully within the boundary. The use of renewable energies, energy efficiency procedures, and also less use of energy here outlined as CHB are the sub components which are located in the overlapping area of the three domains of sustainability. Based on this diagram a general visual representation of policy4A.2 analysis is demonstrated using a Probabilistic Network in Fig. 6: The above general model consists of five modules which will be further detailed in the context of this work. The decision node or policy node is '30% Carbon Footprint (CF) reduction by 2025'. This is one of the nodes that do not include a mathematical expression while it only demonstrates an information passing process, rather than actual calculability.

The above decision node includes: The linearity of the output is due to the input being the contents of a table of data rather than a non linear equation. The probability distribution of the cost for not investing is as shown below: If Projected_cost>Total_Costs, Then True, Else False

The above influence diagram and its modules and sub diagrams are to be created after this initial generalized model. But the reasons for choosing the influence diagrams for this type of policy modeling are due to its advantages. Influence diagrams have many advantages over decision trees. In decision trees adding nodes and generally expanding the model with more details requires to add nodes as child within the tree and that adds nodes to the whole diagram exponentially. Another issue with that is articulating conditional independence relationships within a decision tree would be easier said than done. Some levels of expressions within influence diagrams cannot be developed using decision trees [START_REF] Howard | Influence diagrams[END_REF].

A key question here would be that of how to evaluate the ID created. There are different naming used for different types of ID. A regular ID is the one with decision nodes in order, one after the other within a specified path. The not-forgetting ID is the one that information for previous decision nodes is also available for other decision nodes. A well-formed influence diagram or abbreviated as WFID is an influence diagram with well defined mathematical and probability expressions within itself. There are different ways in order to evaluate an influence diagram. The evaluation technique for WFID would be the reversal and removal operations in order to reach the starting node. For regular and other forms of influence diagrams there is the chance for conversion of the ID into a Bayesian network, although these two are quite the same. Converting an ID into a Bayesian Network is a straightforward task, as both have many similarities. All the nodes within the ID should be converted into chance nodes. The reason for that is in a Bayesian network the nodes are representing a probability distribution and no exact node should be defined. For conversion of the decision nodes into probability nodes or chance nodes, the newly defined chance node should have an even distribution:

∀di ∈ dom(di),P(di |πdi) = 1 /| dom(di)|
Where di ∈ D corresponds to decisions to be made It is once simulated for a probability less than 50% (can be tentatively named true chance), and once for a probability over 50% (which alternatively can be named as the false chance). Then it can be simulated as Bayesian network and evaluated using the Bayesian algorithms. This technique is called Coppers reduction which was initially introduced by Cooper in 1988. There are other algorithms and techniques which will enable us to evaluate ID e.g. Shachter and Peot's Algorithm, Zhang's Algorithm, Xiang and Ye algorithm and etc [START_REF] Crowley | Evaluating Influence Diagrams: Where we've been and where we're going[END_REF]. The completed ID for policy 4A.2 can be evaluated by means of an extensive comparison between the result of the ID and a scientific project management approach. The evaluation methods discussed above can also be utilized for that purpose.

Conclusion and Further Work

It appears that the concept of Influence Diagrams is compatible with the nature of decision making challenges paused by the London Plan and in particular the financial viability of the use of Renewable Energy Systems to partly meet LP's targets. Completion of the above ID would be a proper means for evaluating the financial domain of the London Plan 4A.x policy. Apart from completing the ID, evaluation of the ID is essential. Policy making for tackling climate change, energy crisis and a probable catastrophe has been urged by many governments. Some have started it already using various technical and scientific research studies. However, some policies are still at political level and there is no proper research work to allow for the validation of such plans. According to the United Nations MDG indicators datasets, China and United States are the most CO2 polluters. Political and cultural similarities between the UK and the USA confirm that the next step for the use this climate change policy modelling approach would be the US cities with the highest CO2 emission [20]. Similar policies can be applied to them. Finally, use of other renewable energy systems would not cause any amendments to the general structure of the research using Influence Diagrams.

Fig. 4 :

 4 Fig.4: Artificial Intelligence Classification[START_REF] Hopgood | Intelligent Systems for Engineers and Scientists[END_REF]

Fig. 5 .

 5 Fig. 5. Systems Map for London Plan Policy 4A.2

Fig. 6 .

 6 Fig. 6. High level Influence Diagram for London Plan Policy 4.A.2

Fig. 7 .

 7 Fig. 7. Decision node: 30% CF Reduction by 2025 considering the 1990's base

Fig. 8 .

 8 Fig. 8 . Cost of Carbon Footprint by 2030, a tentative output illustrating different possible degrees of membership and the CF cost.

Table 1 .

 1 Energy Efficiency and Renewables savings in the London Plan

		Savings to date from 113 energy	Savings to date -scaled for 350
		statements		developments
		tonnes CO 2 /y	%	tonnes CO 2 /y
	EE	111,492	21.3	345,329
	RE	24,039	5.8	74,457
	Total	135,528	25.8	419,777


Ping-Ho Chen

Kuang-Yow Lian
email: kylian@ntut.edu.tw

Optimal Fuzzy Controller Mapped from LQR under Power and Torque Constraints

Keywords: LQR, FLC, Fuzzy, Sugeno, Optimal, Class stacking

Dealing with a LQR controller surface subject to power and torque constraints, is an issue of nonlinear problem that is difficult to implement. This paper employs a fuzzy controller surface to replace the LQR surface subject to power and torque constraints by using class stacking, least square and Sugenotype fuzzy inference mode. Through this type of transformation, called "Optimal fuzzy controller mapped from LQR", control of the system remains optimal.

1

Linear Quadratic Regulator (LQR)

The development of optimal control for a MIMO system [START_REF] Frank | Applied Optimal Control & Estimation[END_REF][2] is mature nowadays. Approach of linear quadratic regulator (LQR) offers Kalman gain for state feedback to regulate a MIMO system in state-space form. But subject to power and torque limitation of motor drive, implementation of system becomes difficult using LQR approach. This section introduces some basis for LQR approach for further improvement.

LQR Approach

For a MIMO

- Bu X BK A X + = ) ( & (3)
At steady state it satisfies

com - Bu X BK A O + = ) ( (4)
In case o u = com in Eq. (4), the system becomes a regulator problem further called "linear quadratic regulator (LQR)", if we select a gain K , called Kalman gain such

that KX u - =
to minimize the cost J in quadratic form as

)] - ( ) ( 2 1 [ ) ( ) ( 2 1 T T 0 T T dt T T J T X Bu AX Ru u QX X X S X & + + + + = ∫ λ (5)
where Q: dynamic state weighting, R: control weighting : S terminal state weighting and λ: Lagrange multiplier. After several algebra operations to minimize the cost, we have the following Hamiltonian matrix presenting as the transition matrix in Eq. ( 6):

⎥ ⎦ ⎤ ⎢ ⎣ ⎡ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ = ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎣ ⎡ λ X A B BR A λ X - Q - - T T -1 & & ( 6
) ) ( ) ( ) ( t t t X S λ = (7)
) ( ) ( ) ( T T T X S λ = ( 8
) ) ( ) ( - ) ( ) ( - ) ( - ) ( T -1 T -1 t t t t t t X K X S B R λ B R U = = = (9)
The steady gain is

-1 11 21 -1 T -1 ) ( ) ( lim ) ( lim ψ Bψ R S B R K K = = = ∞ → ∞ → t t t t (10)
where 21 ψ and 11 ψ are the first-column block matrix of Ψ which is the modal form of the Hamiltonian matrix in Eq. ( 6).

Constraints of LQR Control

While implementing the optimal control u aforementioned, we always meet difficulty when subject to constraints of torque limit l u and power limit l p . The torque constraint, caused by stall current in motor application, is expressed by

s i l l k i u u u = ≤ , ( 11
)
where i k : current constant of motor. s i : stall current of motor. As usual the control u is associated with states of position and rate. Mechanical power m p , constrained by maximum power l p of motor, is represented by

l m p p ≤ = θ u & T ( 12
)
where θ & is rate-related state vector.

Class Stacking

From Eq. ( 2), the LQR control is

t 2 1 ] [ w u u u - L = = KX u where n w ji R k × ∈ = ] [ K and ∑ = = n i i ji j x k u 1 - for w j , 1,L =
. Since each component j u of u will be solved one by one by the same approach of class stacking,

Fuzzy inference system (FIS)

A fuzzy logic controller (FLC) [START_REF] Zadeh | The concept of a linguistic variable and its application to approximate reasoning[END_REF][5] [START_REF] Mamdani | Advances in the linguistic synthesis of fuzzy controllers[END_REF] is configured by processing a sequence of fuzzification, inference and defuzzification as follows. The fuzzy inference system (FIS) to be selected is Sugeno model [START_REF] Sugeno | Fuzzy measures and fuzzy integrals: a survey[END_REF] given in Fig. 2 with 1 υ and 2 υ as input fuzzy variables and the control u as known output variable which is bounded, constrained and expressed by plane s1:

Fuzzification

Input fuzzy variables are 1

Inference

Output variable is the control u that is composed of plane s1 and a hyperboloid.

The plane s1 is Λ is obtained by

) ( ) ( T -1 T T Z Y Y Y Λ = (21) with ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ = ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ⎡ = D C B A D 1 D 2 C 1 C 2 B 1 B 2 A 1 A 2 u u u u Z Y , 1 1 1 1 υ υ υ υ υ υ υ υ (22)
The selected Λis able to minimize error YΛ -Z .

After finishing the definition of linear planes, we are able to develop the inference that includes ten rules. The former half governs the surface s2 and the latter half governs the surface s1. Rules are written as follows.

1. If ( 1 υ is S) and ( 2 υ is B) then ( u is s2); 2. If ( 1 υ is SM) and ( 2 υ is MB+) then ( u is s2) 3. If ( 1 υ is M) and ( 2 υ is M+) then ( u is s2); 4. If ( 1 υ is MB) and ( 2 υ is SM+) then ( u is s2) 5. If ( 1 υ is B) and ( 2 υ is S+) then ( u is s2); 6. If ( 1 υ is S) and ( 2 υ is not B) then ( u is s1) 7. If ( 1 υ is SM) and ( 2 υ is not MB+) then ( u is s1); 8. If ( 1 υ is M)

. Defuzzification

A rule surface after defuzzification is generated as shown in Fig. 3(b) which is observed to be an approximation of Fig. 3.(a).

Conclusion

LQR control surface in hyperspace, subject to power and torque constraints, is difficult to implement. This paper employs class stacking to transform LQR surface in hyper space into a surface in 3-D and further applies a Sugeno-type fuzy logic controller to plot the 3-D nonlinear surface. This approach makes LQR feasible and visible.

(

  without loss of generality. Obviously for j th component j u of u , plane of n-D. All the states are sorted out into two classes. One class, having m states, is position-related and piled up by class, having n-m states, is rate-related and piled up by 13)~(15) thus reduce the system from n-D into 3-D with variables 1 control u but upper bounded by l u . Eq. (17) indicates a nonlinear constraint. The process, called "class stacking", thus make it feasible to plot a hyper plane of n-D in terms of 3-D. An example using computer-aided plot is given as follows. A LQR control u , with number of total states n=5 and number of position-related states m=3, Inequality constraint). By using computer-aided plot [3], Fig. 1(a)-(d) shows a sequence of transformations from LQR surface to cases of constraints are obtained. Fig. 1(a) shows the interception of a LQR surface, a torque-constrained surface and a hyperboloid due to power constraint. Fig. 1(b) shows LQR surface bounded by all constraints.

Fig. 1 .

 1 Fig. 1. Computer-aided plot (a) Interception of LQR plane, torque constraint plane and power constraint hyperboloid (b) LQR surface bounded by all constraints.

  Although the approach of class stacking reduces the hyper plane from the LQR controller of n-D into a bounded and constrained plane in 3-D, the 3-D controller surface shown in Fig. 1.(b) is still nonlinear. Therefore we may employ a fuzzy-logiccontroller (FLC) to implement this 3-D controller surface shown in Fig. 1.(b).

υ and 2 υ . 1 υ

 1 Fig. 2.

Fig. 2 .

 2 Fig. 2. Fuzzy inference system (FIS) using Sugeno model

.

  The hyperboloid ABCD in Fig.1(b) is approximately replaced by another plane in Fig.3.(a), named s2 and obtained by least square with s2:

Fig. 3 .

 3 Fig. 3.(a) Hyperboloid is approximately replaced by plane s2 (b) Rule surface after fuzzy defuzzification.

  system in state-space form

	Gaussian) problems. When	u = dist	o	and	u ≠ com	o	,the system becomes a tracking
	problem expressed by
												com
								X &	=	AX	+	Bu	(1)
	where	X	∈	R	n	is the state vector. The control	u	∈	R	w	includes state-feedback torque in
	motor application, command for desired state output and disturbance that the system
	might encounter. Control u is thus expressed by
						u	=	-KX	+	com u	+	dist u	(2)
	where K : Kalman gain. dist u : disturbance torque (colored noise). com u : command
	for desired state output. Eq. (2) is a general form of input in LQG (Linear Quadratic