Yulin Zhang

Yang Xu
email: xuyang@uestc.edu.cn

Haixiao Hu

Cooperative Decision Algorithm for Time Critical Assignment without Explicit Communication

Keywords: Multi-agent systems, task assignment, time critical

Decentralized time critical assignment is popular in the domains of military and emergency response. Considering the large communication delay, unpredictable environment and constraints, to make a rational decision in time critical context toward their common goal, agents have to jointly allocate multiple tasks based on their current status. This process can be built as a Markov Decision Model when they can share their status freely. However, it is computationally infeasible when explicit communication becomes a severe bottleneck. In this paper, by modelling the decision process as a Partially Observable Markov Decision Process for each agent and building heuristic approaches to estimate the state and reduce the action space, we apply greedy policies in our heuristic algorithms to quickly respond to time critical requirement.

Introduction

With the progress of research in distributed artificial intelligence, cooperative agents have been applied in dynamic and complex application domains such as disaster response [START_REF] Scerri | Flood disaster mitigation: A real-world challenge problem for multi-agent unmanned surface vehicles[END_REF], military [START_REF] Day | Multi-agent task negotiation among UAVs to defend against swarm attacks[END_REF] and business organizations [START_REF] Luo | A distributed algorithm for constrained multirobot task assignment for grouped tasks[END_REF] in recent years. In those domains, task allocation becomes a basic problem that needs to be solved efficiently [START_REF] Luo | A distributed algorithm for constrained multirobot task assignment for grouped tasks[END_REF]. The basic task allocation problem is organized as follows: Given a set of agents and a set of tasks, where each agent obtains some payoff for each task, find a one-to-one assignment of agents to tasks to maximize the total payoff of all agents. This problem can be solved (near) optimally in polynomial time by centralized algorithms [START_REF] Kuhn | The hungarian method for the assignment problem[END_REF][START_REF] Burkard | Assignment problems[END_REF], and decentralized algorithms [START_REF] Luo | A distributed algorithm for constrained multirobot task assignment for grouped tasks[END_REF][START_REF] Bertsekas | The auction algorithm: A distributed relaxation method for the assignment problem[END_REF]. In all of these works, it is assumed that communication is available and the environment is static. However, in real domain, the communication is not available and the environment including the status of the agents and the tasks may change. For example, when a group of robots are deployed to carry out time critical tasks for Urban Search and Rescue (USAR), communication infrastructure is most likely destroyed in disaster and communication may be restricted for each rescue robot. Another example, when missiles are launched to attack hostile targets, the missiles have to keep silent due to security reasons. In both of these domains, the environment is dynamically changing and necessary cooperative replanning or reallocation is required. Therefore, we summarize those scenarios as a group of time critical agents that they have to cooperate to optimize their common goals without explicit communication.

More specifically, time critical task assignment without explicit communication that we study can be stated as follows: Given a set of agents and a set of tasks, firstly, the capability of each agent is determined by the dynamically changing environment; Secondly, the state of task and agent may change unexpectedly; Thirdly, no explicit communication is allowed; In addition, a payoff can be obtained from a completed task. We need to find the assignment of agents to tasks within critical response time, such that the sum of payoff from all agents is maximized. To improve team performance, team states have to be well maintained from observations. This can be modelled as a Decentralized Partially Observable Markov Decision Process, which is NEXP-complete [START_REF] Bernstein | The complexity of decentralized control of Markov decision processes[END_REF]. And current approaches such as [START_REF] Nair | Networked distributed pomdps: A synthesis of distributed constraint optimization and pomdps[END_REF][START_REF] Zhang | Coordinated multi-agent reinforcement learning in networked distributed pomdps[END_REF], try to solve this problem with the help of communication, which is not allowed in our problem.

In this paper, several efforts were made to solve the task assignment problem in time critical context. Firstly, since the capabilities of all agents and the statuses of their tasks are changing with Markovian property, the task assignment process can be built as a Markov Decision Model. However, when communication is not available or unacceptable in time critical context, the team state cannot be accurately observed. In this paper, we model the decision process as a Partially Observable Markov Decision Process for each agent. Secondly, the objective function in our task assignment is a step function, which makes it difficult for the agent to decide whether to jointly assign an uncompleted task to a single agent. We propose an exponential reward function to fit for it. Finally, without communication, there is only a small amount of observation, which makes the problem an open-loop control and toward an NEXP-complete problem [START_REF] Papadimitriou | The complexity of Markov decision processes[END_REF]. Thus, heuristic approaches are built according to the following observations: 1. Task-related agents may get the same observation with a high probability. 2. When nothing can be referred from others, it's better to assume that their choice stays the same. With these observations, the search space can be dramatically reduced and a greedy policy is applied to help operate in time critical context.

Problem Statement

In this section, we provide a time critical assignment problem without explicit communication for cooperative teams.

Task Assignment Problem

We focus on a cooperative time critical agent team A = {a 1 , a 2 , ..., a i , ...}, such as a group of robots or UAVs, coordinating to carry out a set of tasks T = {t 1 , t 2 , ..., t j , ...}. Each task consists of two parts, t j =< χ j , Φ j >. χ j is the payoff to complete t j , indicating its importance. Φ j describes t j 's workload, which represents the total amount of work required to complete t j . Before the task is executed, there is a preparation stage for the agent to prepare for the execution or heading for the target location. During this preparation stage, the potential contribution of the agents to their tasks may change, since their states dynamically change. To obtain the highest payoff with limited agents, we should dynamically assign all the tasks during this stage.

At time t, the potential contribution of a i to t j is described as c(a i , t j ) t , which is dynamically updated as the situation between a i and t j changes, according to the domain knowledge. For example, the capability of the rescue robot is related to its remaining energy, holding resources and even distance to its target location. When the robot is of low power, its capability to carry out a long-distance task is very small. If the resource in the robot is going to run out, it is less capable to carry out a resource-consuming mission. After agent a i is assigned with task τ ai (t), all agents assigned with task t j is summarized as Γ j (t). Thus, if the total contribution ψ j (t) of all assigned agents for task t j meets its workload, that is,

ψ j (t) = ∑ ai∈Γj (t) c(a i , t j ) t Φ j t j can be completed and Completed(t j ) = 1; Otherwise, Completed(t j ) = 0.
During the preparation stage, we need to find the best strategy π * as following to get the highest payoff.

π * (t) = argmax Jt τ (A,T,t) ∑ tj ϵT χ j • Completed(t j ) (1)
where Jt τ (A, T, t) is the task allocated for all agents in this team, Jt τ (A, T, t) = ∪ ai∈A τ ai (t).

Time Critical Assignment without Explicit Communication

In a time critical system, a tiny delay may lead to disasters. For example, when a UAV in fast mode finds an obstacle in front, the decision about how to change its orbit should be made immediately, or it will be destroyed before it changes its orbit. Therefore, the task assignment should be carried out in a time critical context and no deliberation model should be allowed. In addition, when the time critical system is performing a task far away from its base, communication is extremely unreliable and of long delay. To operate in a time critical context, no explicit communication could be relied on. Without explicit communication, due to the limited sensing capabilities, the agents cannot fully observe the team state but can only get partial observations. Therefore, the agents have to reason from its partial observations to make rational decisions.

Cooperative Task Assignment Algorithm

In order to obtain the highest payoff, each agent has to dynamically update its task according to the changing team state, including the capabilities of all agents and the status of each task. Since the task of each agent is only determined by the latest team state rather than the historical ones, the task assignment process can be built as a Markov Decision Model when they can share their status freely. Unfortunately, in time critical context, the agents cannot share the state through unreliable communication with large delay. Each agent can only estimate the team state from the ambient data of its own observation. Due to the limitation of its observation range, the observation obtained by the agent is incomplete and team state is uncertain. Therefore, in this paper, a Partially Observable Markov Decision Process model is built for each agent to reason from its partial observations. Since small amount of observation makes the POMDP problem even more difficult, heuristic approaches are proposed to solve this problem.

Partially Observable Markov Decision Process Model for Cooperative Task Assignment

The basic decision model for the task assignment process of each agent can be written as a POMDP < S, Jt act, T, Θ, O, U >. S is the state space and its value at time t is defined as s(t), Jt act is the action space of all agents, T : S ×A → S, is the transition function that describes the resulting state s(t + 1) ∈ S when executing Jt τ (A, T, t) ∈ Jt act in s(t). Θ and O are the observations and their observation function respectively to indicate the possible states under a given observation. R : S → R defines the reward of being in a specific state. This model can be applied to each agent. In this case, the team state s(t) is modelled as the capabilities of all agents and the assigned agents for each task at time t -1:

S(t) = ∪ tj ∈T (Γ j (t -1) ∪ ai∈A c(a i , t j ) t )
After assigning Jt τ (A, T, t) at state s(t), the capability of each agent is transited according to the domain knowledge: T : c(a i , t j ) t ×τ ai (t) → c(a i , t j ) t+1 . For example, when a rescue robot changes its task, it will change its direction to head for its new target location and improve its capability to accomplish this task in time.

However, each agent cannot fully observe the whole team state but can only get partial observations θ(t) ∈ Θ, such as the status of the task, the state of itself and other agents. For example, when the rescue robot is approaching its target location, it can sense the topography of the target area and its remaining energy. From the observation, the belief state b(t) could be derived according to the observation function: O : Θ → B. For each state s(t) ∈ b(t), we know the perceptual distribution P r(s(t)|θ(t)), which describes the likelihood of the team being in the state of s(t) under observation θ(t). For example, when a robot observes that the traffic to the target location becomes extremely tougher, which will consumes more energy to reach the target location and indicates a capability drop of this agent, others with the same task may also observe the same capability drop with a high probability.

To find the optimal assignment, the key is to maximize the objective function [START_REF] Scerri | Flood disaster mitigation: A real-world challenge problem for multi-agent unmanned surface vehicles[END_REF]. However, the objective function is a step function, which generates reward only when the task is able to be completed and makes it difficult to decide whether to assign an uncompleted task to a single agent. Hence we propose an exponential reward function to fit for the objective function. And the exponential reward function should satisfy the following two properties: Firstly, for each task t j , with all the potential agents and their contributions, the potential progress ψ j (t) of t j is the sum of all contributions from its assigned agents.

ψ j (t) = ∑ ai∈Γj (t-1)
c(a i , t j ) t If ψ j (t) equals to zero, which means that there is no potential progress for task t j and the team obtains no reward from task t j . Secondly, when the task receives enough progress to be completed, the team gets the whole priority of the task as its reward. Considering these two properties of the objective function, the exponential reward function is built as follows:

R(ψ j (t)) = χ j e Φj -1 (e ψj (t) -1)
And the utility of s(t) is the total reward from all tasks.

U (S(t)) = ∑ tj ∈T R(ψ j (t))
In this model, the agents coordinate to find the best strategy π * to maximize the expected utility.

π * (t) = arg max Jt τ (A,T,t) ∑ s(t+1)∈b(t+1)

P r(s(t + 1)|θ(t + 1)) × U (s(t + 1))

Without communication, each agent can only get a small amount of observation, which makes this problem an open-loop control and toward an NEXP-complete problem, according to the study of [START_REF] Papadimitriou | The complexity of Markov decision processes[END_REF].

Heuristic Approach

In this section, we build heuristic approaches according to the following efforts. Firstly, we assume that an optimal assignment is initialized for each agent offline. Secondly, since agents act in the same environment, task-related agents may obtain the same observation with a high probability. Hence the belief of all taskrelated agents may change and their tasks may be reassigned. Finally, when we cannot infer anything about others, since the state is Markovian, we can assume that their choice stays the same to keep synchronization of all agents. Thus, the search space can be dramatically reduced and a greedy policy can be applied for each agent in a time critical context.

State estimation.

In a decentralized task assignment domain, several agents may be assigned to the same task. When there is an unexpected accident happening to the task, all agents assigned to the task tend to get the same observation. For example, when a group of robots are assigned to carry out the same task t, if robot α observes that the condition of the task is more difficult than it expected, its capability c(α, t) t will drop and other task-related agents may also observe the same disturbance, which drives them to reassign their tasks. Thus, the state can be estimated with all the possible observations from all task-related agents.

s(t) = s(t -1) ∪ ai∈Γt(t-1) c(a i , t) t
where t is the task that causes the disturbance and c(a i , t) t is the estimated capability from the observations. Action space reduction. In state estimation, task-related agents in Γ t may get the same observation and their states can be easily estimated. But without communication, we know nothing about other agents, which do not share the same task and their states cannot be derived. Hence there are many possible states for these agents, which make it difficult to reason about their states within critical response time. When we know nothing about other agents and what they know, since the state is Markovian, it is reasonable to assume that their state changes regularly and their choices stay the same. That is, the joint action of the agents in A -Γ t , whose state cannot be derived, is estimated from that of last assignment.

Jt τ (A -Γ t , T, t) = Jt τ (A -Γ t , T, t -1)
Thus, only the joint action Jt τ (Γ t , T, t) of the task related agents in Γ t has to be reallocated, which dramatically reduces the action space. With these two parts, the joint action of the team is approximated:

Jt τ (A, T, t) = Jt τ (Γ t , T, t) ∪ Jt τ (A -Γ t , T, t)
To react in a time critical context, no deliberation model could be used. Therefore, we use greedy policy to find a good joint action Jt τ (Γ t , T, t) for the task-related agents.

Each agent in the team runs algorithm 1. In this algorithm, when the agent gets an observation (line 2), it triggers a task assignment process. For the agents that have different tasks (line 4), it assumes that their tasks keep the same (line 5). According to the domain knowledge, their contributions to the tasks can be updated (line 7). Thus the total potential contribution each task receives could be summarized (line 8). For the task-related agents (line 10), the contribution to each task can be updated from the observation (line 11). To obtain the highest priority and respond within critical time, a greedy policy is conducted to find the task with max reward for each agent. Considering the task may be over completed, the reward has to be calculated in two cases (line 14-19 and line 20-26). When the task is not over completed, the reward is the additional utility after adding current agent (line 15). If the reward is the maximum one, record the task and its reward (line 16-17). When the task is over completed after adding current agent (line 20-21), the reward is the additional utility of the task before it exceeds its priority (line 21). And the same is done to record the maximum reward (line 22-23). Finally, find the task with maximum reward and update its status (line 28).

Algorithm 1 Cooperative Decision Algorithm for Time Critical Assignment

Experiment Setup and Result

In this section, we present the empirical evaluations of our approach in the cooperative task assignment. In our experiment, a group of agents are deployed from the same base station to carry out several tasks, which takes the agents about 40 time steps to reach the target location. Each task is initialized with different priorities and workloads. For the sake of simplicity, we initiate agents' capabilities for all tasks to 1. The objective is to maximize the priorities of all tasks accomplished by the agent group. Initially, the agents do not have accurate knowledge about all the tasks and the situation may change. When the agent is heading for the target location, its capabilities can be found changed since the situation is different.

In order to investigate the capability of adapting to the dynamics and uncertainties, we compare our approach (dynamic approach) with traditional static assignment policy (static approach), which initially assigns the optimal policy to each agent offline. Since the agents are initialized with identical capabilities, 0-1 Integer Programming could be used to find the initial optimal policy. To explore the performance of our approach, we initialized 20 agents to carry out 20 tasks with different priorities, and compared our approach with static approach and the upper bound. During this experiment, 5 tasks change their resource requirement at different time. Thus, when the agents approach their target locations, they will find their capabilities changed. According to the experimental result in Figure 1, both the total priorities of dynamic approach and static approach drop when there are unexpected disturbances. However, since dynamic approach reasons from the observations and reassigns its task, there is less performance loss in dynamic approach than in static approach.

In order to investigate the scalability, we vary the number of agents and tasks as 20, 50, 100, and 200 in the experiment in Figure 2(a), and also conduct the experiment with different numbers of disturbances as shown in Figure 2(b). The performance is presented as the ratio of the obtained priorities to the upper bound. From Figure 2(a), when the numbers of agents and tasks scale up, given number of disturbances can only influence a smaller part of tasks. Thus, these disturbances have smaller influences on the team performance when the team and tasks scale up. But dynamic approach always performs better than static approach. When the number of disturbances increases as shown in Figure 2(b), more tasks are influenced and the team performance drops. Because of the reassigning capabilities, dynamic approach still outperforms the static one.

Conclusion

In this paper, we presented a multi-agent decision algorithm for time critical assignment without communication. In this context, each agent can only obtain a partial observation of the team state and cannot communicate with each other. To reason with partial observations, we built a POMDP model for each agent and proposed heuristic approaches to estimate the state and reduce action space. Our experiments show that our algorithm is feasible to dynamically adjust to disturbances and can be applied in time critical assignment without communication.

Fig. 1 .

 1 Fig. 1. The experimental results of 20 robots to carry out 20 tasks

Fig. 2 .

 2 Fig. 2. The experimental results of different scales and disturbances


Ping Jiang

Quansheng Dou

Automated Localization and Accurate Segmentation of Optic Disc based on Intensity within a Minimum Enclosing Circle

Keywords: Automated localization, Accurate segmentation, Intensity threshold, Minimum enclosing circle

This paper presents a method for automated localization and accurate segmentation of the optic disc. An intensity threshold is determined and select all the pixels whose intensities are greater than the threshold, by erosion the optic disc can be localized. By dilation and region filling, a minimum enclosing circle which can completely hold the optic disc is determined, within the circle, the vessels are eliminated by replacing the darker vessel pixels with brighter pixels. Define the intensity features of the optic disc boundary , and select the pixels according to the features, then the optic disc may be segmented. The experiment shows that compared to the active contour models, this method is more efficient and accurate on the boundary extraction of the optic disc.

Introduction

The fundus images are used for diagnosis by trained clinicians to check for any abnormalities or any change in the retina. The information about the optic disc can be used to examine the severity of some diseases such as glaucoma. Its detection is prerequisite for the segmentation of other normal and pathological features.

A number of studies have reported on automated localization of optic discs; several studies have also reported on the segmentation of optic discs [START_REF] Walter | Segmentation of color fundus images of the human retina: detection of the optic disc and the vascular tree using morphological techniques[END_REF][START_REF] Lalonde | Fast and robust optic disc detection using pyramidal decomposition and Hausdorff-based template matching[END_REF][START_REF] Li | Automated feature extraction in color retinal images by a model based approach[END_REF][START_REF] Osareh | Comparison of colour spaces for optic disc localization in retinal images[END_REF][START_REF] Lowell | Optic nerve head segmentation[END_REF][START_REF] Wong | Level-set based automatic cup-todisc ratio determination using retinal fundus images in ARGALI[END_REF][START_REF] Abramoff | Automated segmentation of the optic disc from stereo color photographs using physiologically plausible features[END_REF][START_REF] Jaspreet Kaur | Automated Detection of Vascular Abnormalities in Diabetic Retinopathy using Morphological Thresholding[END_REF][START_REF] Sinthanayothin | Automated localisation of the optic disc, fovea and retinal blood vessels from digital colour fundus images[END_REF][START_REF] Hoover | Locating the optic nerve in a retinal image using the fuzzy convergence of the blood vessels[END_REF][START_REF] Park | Locating the optic disc in retinal images[END_REF][START_REF] Herbert | Automated Image Detection of Retinal Pathology[END_REF][13][14]. Walter and Klein [START_REF] Walter | Segmentation of color fundus images of the human retina: detection of the optic disc and the vascular tree using morphological techniques[END_REF] applied the watershed transformation to the gradient image based on the morphological operations. Lalonde et al. [START_REF] Lalonde | Fast and robust optic disc detection using pyramidal decomposition and Hausdorff-based template matching[END_REF] used the Canny edge detector to detect optic disc edge, and by matching the edge map with a circular template the optic disc was segmented. The method of Li and Chutatape [START_REF] Li | Automated feature extraction in color retinal images by a model based approach[END_REF] is based on an active contour model, and it iteratively matching the landmark points on the disc edge. Osareh et al. [START_REF] Osareh | Comparison of colour spaces for optic disc localization in retinal images[END_REF] and Lowell et al. [START_REF] Lowell | Optic nerve head segmentation[END_REF] used the active contour model to find the optimal points based on the external and internal energy of the image. Wong et al. [START_REF] Wong | Level-set based automatic cup-todisc ratio determination using retinal fundus images in ARGALI[END_REF] proposed a method based on the level-set technique, and used ellipse to fit the disc boundary. Abramoff et al. [START_REF] Abramoff | Automated segmentation of the optic disc from stereo color photographs using physiologically plausible features[END_REF] employed a pixel classification method using the feature analysis and nearest neighbor algorithm, then group each pixel to rim, cup or background. Jaspreet Kaur et al. [START_REF] Jaspreet Kaur | Automated Detection of Vascular Abnormalities in Diabetic Retinopathy using Morphological Thresholding[END_REF] and Siddalingaswamy P. C. [START_REF] Sinthanayothin | Automated localisation of the optic disc, fovea and retinal blood vessels from digital colour fundus images[END_REF] employed iterative thresholding method followed by connected component analysis to automatically localize the approximate center of the optic disc, then the geometric model based implicit active contour is employed to obtain accurate optic disc boundary. Hoover and Goldbaum used a ''fuzzy convergence'' algorithm to correctly identify the optic disk location in 89% of 81 images [START_REF] Hoover | Locating the optic nerve in a retinal image using the fuzzy convergence of the blood vessels[END_REF]. The method of Park et al [START_REF] Park | Locating the optic disc in retinal images[END_REF] found the brightest pixels by employing the repeated thresholding technique, then used the roundness of the object to detect optic disk features, and then localized the optic disk by using the Hough transform.

Yet because of fuzzy boundaries, inconsistent image contrast or missed edge features, it is difficult to accurately localize and segment the optic disc. Based on the intensity and shape features, this paper proposed a method for the automated localization and segmentation of the optic disc, without human intervention, the optic disc may be segmented accurately and efficiently by the contour extraction, the contour can then be used to help the doctor for further analysis.

Appearance intensity threshold determination for optic disc localization

The optic disc is often the brightest object, ranging from white to yellow, is circular and of reasonably consistent size from patient to patient [START_REF] Herbert | Automated Image Detection of Retinal Pathology[END_REF]. So the pixel whose intensity is brighter than the threshold is marked as candidate optic disc pixel. For the retinal images whose sizes are 640×480, the optic discs have more than 1000 but less than 10000 pixels. Clustered those pixels whose intensities are equal when divided by a certain histogram bin width which in this paper is 10, and sort the clusters according to the intensity. Start from the cluster with maximum intensity and check the number of pixels in it, if the number is less than 1000, then go to the next cluster in descending order of intensity until the number of pixels is between 1000 and 10000, and the intensity range of this cluster is determined as the threshold T . Those pixels whose intensities are brighter than the threshold are marked as candidate optic disc with a green color (0xff00ff00) as shown in Figure1.

Figure1. Retinal image with candidate optic disc marked with green color

As shown in Figure1, there are some regions marked with green color which are not real optic disc, to eliminate the false candidate pixels, morphological erosion is applied on them. Find all the start positions of the candidate regions, structure element structElem begins from these positions and tests all the candidate regions so as to find a region which can include it. Figure2 shows an instance of the localization result.

Figure2. The result of automatic localization of the optic disc with a rectangle

Dilate the optic disc region, determine the minimum enclosing circle

As shown in Figure2, the green pixels don't occupy the whole disc region, so further steps are needed to gather more optic disc pixels.

Morphological dilation

By region grouping, the optic disc pixels are collected together and marked in binary matrix M . The matrix is the same size as the retinal image, if the pixel on the image is colored by green, then the corresponding cell in the matrix is 1, otherwise 0. To find more optic disc pixels, double dilation was applied on M , The dilation result was shown in Figure3 with blue color, take all the x and y coordinates of dilated matrix M , and compute the average

x and y as the coordinate of the center point ( , ) P x y of the blue region.

Figure3. Dilation result of Figure2.

Determination of minimum enclosing circle

Based on the normal size of optic disk, centered at P , a circle which can totally enclose the optic disc region can be determined. Analyze the intensity distribution of pixels in the circle, and define the intensity features of the pixels on the edge, and mark those pixels having edge features, then the segmentation is done. Two edge features are defined initially: the first is that its intensity is darker than the optic disc localization threshold T ; the second is at least one of its 8 neighbors belongs to the internal optic disc, i.e. the intensity is brighter than the threshold T . According to the two features, some pixels are selected which are used to compute their distances to the center point pCenter , and get the average distance avgD , set the radius of minimum enclosing circle to be D is calculated by equation ( 1). The determined edge result is shown in Figure4.

Figure4. The optic disc edge results with the third feature added.

Experiment and Conclusion

The automated localization and segmentation of the optic disc is implemented by JAVA, and it can accurately and efficiently get the contour of the optic disc, yet for the active contour model, to get the optic disc contour, an initial contour need to be drawn first, and because of the noise, vessel intervention etc., normally they cannot get to the accurate boundary, so compared to these methods, our algorithm is more accurate and efficient.

We test our method on DRIVE, DIARETDB1 and the infant retinal images, for segmentation, the method has good performance on the infant images whose optic discs have relatively weak vessel intervention, for those images with strong inside vessels, the method cannot segment all the optic disc, normally a large part the disc, so it needs further improvement which is our future work.



  , where d is the incremental constant used to make the circle include all the optic disc pixels and exclude the false optic Automatic localization and segmentation of optic disc 5 disc pixels. So the third feature of the edge pixel is that its distance D to pCenter should satisfy || D avgD d  .

Acknowledgements

This work was partially supported by National Natural Science Foundation of China, under grant No. 61272244, 61175023,61175053 and 61272430, Natural Science Foundation of Shandong Province, under grant No. ZR2013FL022, the Young Foundation, under grant No. 2011QN083.

AUTHORS' INFORMATION


Chaojun Wang
email: wangchaojun009@126.com

Zhixin Li

Canlong Zhang
email: clzhang@gxnu.edu.cn

A Multi-instance Multi-label Learning Framework of Image Retrieval

Keywords: multi-instance multi-Label learning, image retrieval, multi-points diverse density algorithm, image retrieval

Because multi-instance and multi-label learning can effectively deal with the problem of ambiguity when processing images. A multi-instance and multi-label learning method based on Content Based Image Retrieve ( CBIR) is proposed in this paper, and the image processing stage we use in image retrieval process is multi-instance and multi-label. We correspond the instances with category labels by using a package which contains the color and texture features of the image area. According to the user to select an image to generate positive sample packs and anti-packages, using multi-instance learning algorithms to learn, using the image retrieval and relevance feedback, the experimental results show that the algorithm is better than the other three algorithms to retrieve results and its retrieval efficiency is higher. According to the user to select an image to generate positive sample packs and anti-packages, using multiinstance learning algorithms to learn, using the image retrieval and relevance feedback. Compared with several algorithms, the experimental results show that the performance of our algorithm is better and its retrieval efficiency is higher.

Introduction

We live in an era in which multimedia information grows at a speed of geometric index and network technologies are popularized quickly .It is more and more difficult that the images we need are found quickly and effectively from the pictures about the ocean. We can know that text-based image retrieval and content-based image retrieval doesn't solve well the content of the image of the ambiguity through the above introduction, and a complete image includes a lot of small parts. If only a single label does not adequately describe the image, but cannot accurately expressed the users are interested in the parts, which can lead that the part of the image that we are own interested in is not retrieved .The content of this article is to divide the image into several parts, take partial features from several pieces of information as the package of the image generated by examples, and contact the examples of the package with the category of the mage. This paper presents a method of the image retrieval based on multi-instance multi-label study. The method takes use of local color of images and texture features to generate multi-sample packages, and give the appropriate subset of classes for objects, and it is no longer the only category tag. We can get the complex high-level concepts derived by the subsets of these categories marked, so that we improve the probability of recall and precision in the process of image retrieval.

Multi-instance multi-label learning

In multi-instance learning, every training package is made of more than one sample, and the sample has no concept, but each training package has a concept marked.

When we apply this method to image retrieval, we will find that it belongs to a variety of categories for images. If we only use a single category marked words, we cannot fully describe the image, and it will be a problem of semantic gap. Considering the image content with multiple semantics, we should learn that it is actually a map from the sample set to the category tag set on. Considering the image content with multiple semantics, we should learn that it is actually a map from the sample set to the category tag set on. Multi-instance multi-label [START_REF] Zha | Joint multi-label multi-instance learning for image classification[C]//Computer Vision and Pattern Recognition[END_REF] defined as follows, χ represent example of space γ represent the category space, and represent by mathematical symbols follows: Multi-instance multi-label learning: given a set of data {

(X 1 , Y 1 ), (X 2 , Y 2 ),…,(X m , Y m )}.Objective to learn f: χ γ → . 2 2 χ ⊆ i X Is a set of sample { χ χ χ 1 2 , ,..., i i i n i }, χ χ ∈ ij (j=1,2,…,n i ) and χ ∈ Y i
is a s appropriate categories tag of X i .{y i1 ,y i2 ,…,y ili } y ik ∈ γ(k= ,2,…,l i ). n i is the number of samples among X i , l i is the number of tags among Y i .

et of 1

There is learning frameworks named traditional supervised learning (single example, single tag) and multi-label learning (single example, multiple tags) and multi-instance learning (multiple example, single tag) besides the above one of sample learning frameworks in this paper. From the above, it is not difficult to find that MIML framework can be converted into other three examples in some special cases. We know that only by using the appropriate expression can capture the important information of image. However ,dealing with the ambiguity of the image by using MIML can better deal with semantic gap problems .MIML also helps to study the high-level semantic concepts, such as elephants, lions on the picture, trees, and some grass. Through the concept of markers we can infer that this is a picture of a description of Africa. This example is the use of MIML and then study the underlying logic of high-level concepts .In order to play the ability of MIML framework, there are generally two kinds of methods, respectively, for example as Bridges and multilabel learning as a bridge.. This paper is to convert Multi-instance multi-label learning algorithm MIMLBOOST [START_REF] Zhou | Multi-instance multi-label learning with application to scene classi-fication[END_REF][3] into Multi-instance learning.

3

The structure of the package and image retrieval

3.1

The structure of the sample package process Image of the sample bag choosing what kind of method determines whether the image semantic information is complete, this will directly affect the retrieval effect. But Maron and Yang's Multiple-Instance learning by maximizing diverse density [START_REF] Maron | Multiple-Instance Learning for Natural Scene Classification[END_REF] use the dividing method of fixed generated package. This method will cause originally belong to the same concept divide into different packages or belong to different concepts merged into one package. Sample package need to define in advance. These will only contain color information. Colors of different parts may be similar which can lead to low accuracy of image retrieval. The method [START_REF] Yang | Image database retrieval with multiple-instance Learning techniques[END_REF] of Dai H.Blet each image fixed generates 4 sample packages. What kind of effect will be caused by this method, For a simple image that is bound to introduce noise, and compared with the complex image of some of the information will not wholly include .This paper will use a partition method which is a network clustering of Self-organizing feature Map (Self -Organizing feature Map, SOM) [START_REF] Kohonen | Self-organizing maps[END_REF].This article uses SOM clustering image segmentation methods [START_REF] Jiang | SOM based image segmentation[C]. //Rough Sets, Fuzzy Sets, Data Mining, and Granular Computing[END_REF] to divide images, we use the color features as SOM network input to cluster. After go to the discrete points and region merging after some steps such as each region of the image, we respectively calculate the R, G, B, 3 kinds of color grayscale average to the corresponding color features. If we only use color features cannot describe an area, different images in similar parts of color have some differences and different parts of the color may be similar. In order to improve the accuracy of image retrieval, we have to use textures which are some characteristics of the image .Texture is an important means of visual image .Also it describes gray space distribution of image. Texture images can show some information about the arrangement of surface structure, such as roughness, contrast, direction, the line as the degrees, neat and rough degrees, etc. This paper uses co-occurre e matrix [START_REF] Marceau | Evaluation of the grey-level cooccurrence matrix method for land-cover classification using SPOT imagery[END_REF] method to deal with the texture feature extraction. nc If we break the base of decomposition k, then we can get a 3 (k-1) co-occurrence matrix. This article will choose 4 most commonly used statistics out of 10. Also it can reflect the scale of statistics as the texture feature. Finally we add color space R, G, B whose gray value is most characteristic's-occurrence matrix has a feature vector of 4 (k-1)+3 dimension. Regard the regions corresponding eigenvectors as an example. So it can be an image in which all areas of its corresponding sample image corresponding to the package.

Image Retrieval

Figure 1 describes the image retrieval process. There is multi-instance and multi-label learning concept: regard an image as a sample package and the image has a lot of class label. It may be more meaningful in some cases if we understand why an image has some class label, the MIML makes these possible. Using the MIML we can know there is some samples to have such class label. So we can match sample directly with class label images. And we also not only can use image package but also can use the class label between the samples for image retrieval. In the second chapter we introduce four models of machine learning, the special case of the MIML is the other three kinds of machine learning methods. And in this article the process of solving for multiple tag sample learning problems can be converted the problem into a single tag sample learning problems by using MIMLBOOST algorithm. What is the purpose of doing that? We can find the benefits of doing that through the experiment. This article uses the multi-diversity density algorithm (MPDD) [START_REF] Böhm | The k-nearest neighbour join: Turbo charging the KDD process[END_REF] from the image training package to concentrated study some of the concepts for image retrieval. The idea of multi-diversity density algorithm is using density of multiple points of comprehensive information to describe the concept what users are interested in. This algorithm is an improvement compared with the density of the traditional diversity algorithm [START_REF] Mian-Shu | Multi-points diverse density learning algorithm and its application in image retrieval[END_REF]. When using the density of diversity as a measure, we must select the point of maximum density as a target. If the dimension of the feature subspace is higher, we can't find out the point of maximum density by iterating through the whole feature space. The other is the feature space and vector are both distributions, this will lead to exist multiple local minima points. While using the gradient descent method can't avoid fall into local extremum but get a local optimal solution. We know that the use of one example of a sample point is difficult to contain the contents of the sample, this will cause the description of the contents of the sample packets to have certain content, in order to solve these problems it uses the MPDD algorithm in this article. MPDD algorithm calculation process is as follows.

1. The use of colour and texture extracts a features package, and we enter a positive bag and negative bag, to divided the examples and extract the description of the characteristics of the sample. 2. The use of the diversity of classic learning algorithm to find the divers density points, and output the search of the density of the point.

Experimental Result

We do the experiment based on this method, and compare the experiment results with other methods including the Multiple-instance image retrieval [START_REF] Dai Hongbin | A Multi-Instance based approach to image retrieval[END_REF] proposed by Dai Hongbin, Maron and Ratan and Lozano-Perez's multiple-instance learning method known as the diverse density algorithm [START_REF] Chengyang | Image database retrieval with multiple-instance learning techniques[END_REF]. We compare the experiment results using the method of Dai Hongbin, Maron and Ratan's Single blob with neighbours method and the method of Yang and Lozano-Perez.

Experimental Comparison of four methods

HongbinDai and ZhihuaZhou use the method of extracting colour features together with the 24 texture features extracted by 24 filters of Gabor wavelet transformation.

The feature vector has 27 dimensions: 3 colour features R, G, B and 24 texture features. The feature vector of every region is considered as examples.

Maron and Ratan's Single blob with neighbours' method constitutes the matrix by colour blobs after sampling the image. The feature vector with 15 dimensions is calculated from the R, G and B colour features in each colour blob and four blobs connected with it. There are 9 image packages and each contains a 15-dimensional feature vector.

However, Yang and Lozano -Perez's method is to convert the images to grey ones. After dividing into twenty overlapping regions, every region of the image is conducted mirror filtering, transforming and sampling. There are 40 image packages with 64-dimensional feature vector. The results of the experimental are as follows:

This article is for each category of images 10 times experiment, we used precision recall such a way to measure the efficiency of various methods of retrieval. Precision is when system is the process, the system returns results in the target image and all returned images of a ratio. In this experiment, 5 interested in image and not interested in the image of an initial search results are shown in figure2. And in the process of experiment increases 5 interested in image and not interested in image relevance feedback result is shown in figure 3. For the above all kinds of methods for initial retrieval and related feedback average retrieval time are shown in table1. As can be seen from the figure 2 and figure 3, the proposed method takes precedence over these three methods. In this paper, the average time for the initial retrieval and feedback time is prior to the other two methods but longer than another method. The main reason is that MIML learning methods is adopted in this paper. Linking category tag with the sample, the time we use more diversity density algorithm to find the point is less than the former two methods but is more than the density of diversity algorithm. The time this algorithm takes is less than the former two but more than Hongbin Dai multi-instance algorithm. Adopting multi-sample and multi-label algorithm, this article can well describe the image. So through next section in this paper, we learned that the algorithm has better precision and recall radio.

Results in different image library

In order to study the efficiency and results on the different size of image library by this method, we respectively remove 500 and 1000 images from the original image library .I repeated the experiment for the method the article put forward.

From table 2 we can see that, with the increase of image database retrieval time also increases accordingly, but increased more slowly through the experiment we can find the time complexity of this article was close to linear. 5 show the results of image retrieval and the related feedback in different size of image library. We can see the algorithm has smooth retrieval performance with the increase number of experimental image from the experience, and the performance does not appear significant decline. So, the algorithm of this article can be an applicable method to retrieve a large scale image scale image library.

Conclusion

This article puts forward a kind of image retrieval based on multiple sample multilabel learning technology. The benefits of this technology are to reduce the semantic gap and solve the problem of ambiguity of the image. It enlighten us use the underlying semantic concept of high-level semantics. Through the above experiment we find this method is helpful to improve the image retrieval precision rate and recall rate. This article implementation process use the method of MIML as the first step to process image, utilize MIMIL Boost algorithm transforms the problem of multiple sample tags into multiple sample single tag questions, then the image segmentation algorithm based on SOM clustering segment image into multiple regions. We combine each region's color characteristics and co-occurrence matrix as a sample package, then use multi-points diverse algorithm for image retrieval and related feedback. It is suitable for medium-sized image library and compare with proposed method in our experiments, the results show that the algorithm of this article retrieval effect is better than that of Maron, Ratan, Yang and Lozano-Perez etc.This paper also introduced to other characteristics to improve the image of package technology.

In the future research we learn how to put the MIML into the image retrieval studying, reduce the time of recall and precision of image retrieval ,and design a better algorithm to find out the target class image in the shortest time.

Fig. 1 .

 1 Fig. 1. Flow diagram of image retrieval

Fig. 2 .Fig. 3 .

 23 Fig. 2. Initial retrieval results of relevant and 5 irrelevant images

Fig. 4 .Fig. 5 .

 45 Fig. 4. Relevant feedback results of 10 relevant and 10 irrelevant images

  The image dataset this experiment uses has a medium scale. It divides into 10 classes: elephants, flowers, birds, horses, dinosaurs, buildings, and cars, mountains with snow, beaches and plants. The total number of images is up to 1500 and each class has 150 images. We pick 30 images from each class to form the training set as examples. The rest 1200 images form a test set. During each experiment, one class is selected as the target class, and then 5 images which we are care about are selected to generate the positive packet. We select 5 classes randomly from the rest classes and then select an image we are not care about from each of 5 classes randomly to form the negative packet. Now a training set is formed and it includes a positive and negative packet which has 5 images respectively. Here we use more diversity, MPDD density algorithm to learn the training set and retrieve the image in the test set. Then the results are returned by sorting the Euclidean distance of similarity. From the search results, we select the first 5 images which are sorted in target and non-target category, and add the corresponding packet and examples to the training set.

Table 1 .

 1 Performance comparison on the task of automatic image annotation

		Maron and Ratan method	Lozano-Perez andYang method	DaiHongbin instance method	Our method
	The average
	time of Initial	3.2	7.63	1.47	2.36
	retrieves (s)
	Relevance
	feedback time	9.56	8.52	2.92	4.24
	(s)

Table 2 .

 2 Performance comparison on the task of automatic image annotation

		500 images	1000 images	1500 images
	The average time of initial retrieves (s)	1.21	1.56	1.72
	Relevance feedback time (s)	2.86	3.12	3.36
	Figures 4 and figure

Acknowledgement

This work is supported by the National Natural Science Foundation of China (Nos. 61165009, 61262005, 61363035, 61365009), the National Basic Research Program of China (No. 2012CB326403) and the Guangxi Natural Science Foundation (Nos. 2012GXNSFAA053219, 2013GXNSFAA019345) and the "Bagui Scholar" Project Special Funds.

3. Integrated multiple search the information of the density of the point to calculate the density of the point and the distance between the unknown sample and the point. And we calculate the distance between the unknown samples and the density of multiple diversity point through the following formula. Set MPDD algorithm output M (M value unknown) a variety of density points, denoted by = { }.An unknown sample package has N samples, denoted by Then the distance between the unknown sample bag and the density of multiple diversity point is

In this paper we choose more density of points is selected in the sample package. In the process of calculation is to calculate the density of each variety point and test set are package sample the Euclidean distance between, then calculated the distance between the minimum value, then sum and average. Finding out and testing set minimum Euclidean distance are examples of the package contained. We know when described different images, the importance of the various features of feature space is different, if some image has obvious characteristic of image color and texture feature. Then we can give the characteristics of each assigned a weight to represent different characteristics of different degree. The distance between the instances with the target is measured weighted Euclidean distance. Expressed in mathematical formula is as follows:

Image retrieval process is as follows:

─ With the method in the article, we can constitute the image of the sample bag by using the image texture and color. For different image, we can merge sample package with the same texture segmentation area and the color characteristics into sample package with new characteristics by using formula 2 in the article. ─ We can use the K Nearest neighbour (k-Nearest Neighbour , KNN) [START_REF] Maron | A framework for multiple-instance learning[C][END_REF] classification algorithms to classify the training focused sample packages and test focused sample packages. ─ Find the feature vector of the target which is appropriate to users by using formula 1 of the article. Then we can calculate the Euclidean distance between the density point of diversity and all of the samples in the unlabeled sample package. And we can average them as the similarity measure, Order the distance from small to large. After that, present the corresponding image to user. ─ According to the rank of similarity measure, we can choice the interested and uninterested images from the retrieval results and label the corresponding classification of package. After that, put it into training set to get relevant feedback.


Peilong Zeng

Zhixin Li

Canlong Zhang
email: zplhang@126.com

Scene Classification Using Spatial and Color Features

Keywords: Probabilistic Latent Semantic Analysis, scene classification, spatial feature, chain code, KNN-SVM classifier

With the increment of images in modern time, scene classification becomes more significant and harder to be settled. Many models have been proposed to classify scene images such as Latent Semantic Analysis (LSA), Probabilistic Latent Semantic Analysis (PLSA) and Latent Dirichlet Allocation (LDA). In this paper, we propose an improved method, which combines spatial and color features and bases on PLSA model. When calculating and quantizing spatial features, chain code is used in the process of feature extraction. At the same time, color features are extracted in every block region. The PLSA model is applied in the scene classification. Finally, the experiment results between PLSA and other models are compared. The results show that our method is better than many other state-of-the-art scene classification methods.

Introduction

In recent years, image understanding and classification has been frequently researched and widely applied to all kinds of practical systems [START_REF] Shimazaki | Scene classification using color and structure-based features[END_REF]. As an important issue of image classification, scene classification aims to label an image among a set of semantic categories (such as mountains and tall buildings) automatically [START_REF] Bosch | Which is the best way to organize/classify images by content?[END_REF]. For example, images in Fig. 1 can be classified into the category of "coast". It differs from the conventional object classification. Scene classification is an extremely challenging and difficult task because of the ambiguity, rotation, variability and the wide change of illumination and scale conditions of the scene images even if for the same scene category [START_REF] Chen | Evaluating effectiveness of latent dirichlet allocation model for scene classification[END_REF]. The images in Fig. 1 also show that a category may include many coast images with different illumination. What's more, a scene is generally composed of several entities and the entities are often organized in an unpredictable layout [START_REF] Quelhas | Modeling scenes with local descriptors and latent aspects[END_REF]. It is difficult to define a set of properties that would include all its possible manifestations and extract effective common features to classify the images to the same category.

As was mentioned in [START_REF] Bosch | Scene classification using a hybrid generative/discriminative approach[END_REF], there were two basic strategies in the literature about scene classification. The first strategy uses low-level features such as global color, texture histograms and the power spectrum. It is normally used to classify small number of scene categories (city versus open country, etc.) [START_REF] Shimazaki | Scene classification using color and structure-based features[END_REF]. The second one uses an intermediate representation before classifying scenes [START_REF] Li | A bayesian hierarchical model for learning natural scene categories[END_REF][START_REF] Oliva | Modeling the shape of the scene: A holistic representation of the spatial envelope[END_REF] and it has been applied to cases where there are a larger number of scene categories [START_REF] Bosch | Scene classification using a hybrid generative/discriminative approach[END_REF]. Biederman [START_REF] Biederman | Aspects and extension of a theory of human image understanding[END_REF] showed that humans can recognize scenes by considering them as a whole one, without recognizing individual objects. Oliva and Torralba [START_REF] Oliva | Modeling the shape of the scene: A holistic representation of the spatial envelope[END_REF] proposed a low dimensional representation of scenes, based on global properties such as naturalness and openness.

A lot of efforts have been made to solve the problem in greater generality, through design of techniques capable of classifying relatively large number of scene categories in the last few years [START_REF] Bosch | Scene classification via pLSA[END_REF][START_REF] Rasiwasia | Latent Dirichlet Allocation Models for Image Classification[END_REF]. In this paper, we calculate spatial features via chain code method instead of Hough transform in [START_REF] Shimazaki | Scene classification using color and structure-based features[END_REF] first and then combine them with color features, and then use PLSA to classify scene images.

The main contributions of our paper lie in two aspects. One is the application of chain code method to the spatial features' extraction. We use the chain code to describe the shape and spatial information of a scene image. The combination of spatial and color features improve the performance of scene classification. The other is that we use PLSA for learning and SVM and KNN classifiers for classifying. An improved classifier KNN-SVM [START_REF] Zhang | SVM-KNN: Discriminative nearest neighbor classification for visual category recognition[END_REF], the hybrid of KNN and SVM classifier, is used. The experiment results prove the good effect of the classifier.

The next section briefly describes the PLSA model. Then, in Section 3, we describe the classification method by applying PLSA model and KNN-SVM classifier to images. Section 4 describes the spatial and color features used to form the visual vocabulary as well as the feature extraction. The details of experiments and results are provided in section 5. Conclusion is drawn in the 6th section.

PLSA model

Probabilistic Latent Semantic Analysis (PLSA) was proposed by Hofmann in 1999 to solve the problem of ambiguity between words [START_REF] Hofmann | Probabilistic latent semantic indexing[END_REF]. It is a generative model from the statistical literature [START_REF] Hofmann | Unsupervised learning by probabilistic latent semantic analysis[END_REF] which is researched a lot and generates many variations.

In text analysis, it is used to discover topics in a document. Here in scene classification, we treat images as documents and discover categories as topics (such as mountain and road). The model is applied to images by using visual words which is formed by vector quantizing color, shape, texture and SIFT features [START_REF] Bosch | Scene classification using a hybrid generative/discriminative approach[END_REF].

A collection of scene images D=d 1 , ..., d N with words from a visual vocabulary W=w 1 , ..., w V are given. The data in a V×N co-occurrence table are defined as N ij =n (w i , d j ), where n (w i , d j ) denotes how often the word w i occurred in an image d j . In PLSA a latent variable model for co-occurrence data associates an unobserved class variable z∈Z=z 1 , ..., z Z with each observation.

The algorithm of PLSA is as follows: In detail, the KNN selects K nearest neighbors of the new images in the training dataset using Euclidean distance, and then the test image is classified according to the category label which is fit most in the K nearest neighbors; while for SVM classifier, an exponential kernel exp-αd is used, where d is the Euclidean distance between the vectors and α is the learned weight of the training image example [START_REF] Schölkopf | Learning with kernels: Support vector machines, regularization, optimization, and beyond[END_REF]. The perspectives of the KNN-SVM classifier are: The algorithm behaves as KNN classifier and it can easily deal with huge multiclass of the scene images when K is small, while it becomes a SVM model when K is large [START_REF] Zhang | SVM-KNN: Discriminative nearest neighbor classification for visual category recognition[END_REF]. No matter the image dataset is huge or small, the classifier can perform very well.

The algorithm is as follows: ─ Use a crude and simple distance function to find a collection of Ko neighbors.

─ Compute the accurate distance on the Ko sample images and select the K nearest neighbors.

─ Compute (or read from cache if the query repeat) the pairwise distance of the K neighbors and the query.

─ Convert the pairwise distance matrix into a kernel matrix. ─ Apply SVM to the kernel matrix to classify the test image and return the result.

Spatial and color features extraction 4.1 Spatial feature extraction by chain code

Chain code is a compact way to represent the contour of an object for shape coding and recognition. The chain code histogram (CCH) [START_REF] Iivarinen | Shape recognition of irregular objects[END_REF] and minimize sum statistical direction code (MSSDC) [START_REF] Wang | A novel direction chain code-based image retrieval[END_REF] are two methods of all and their advantages are translation and scale invariant. However, they don't take the direction and spatial information into consideration. In the paper, we adopt an advanced method, chain code spatial distribution entropy (CCSDE) [START_REF] Sun | Contour-Shape recognition and retrival based on chain code[END_REF], to calculate the spatial features because it takes full advantage of the statistical feature, the distribution and the relativity of the chain code sequence [START_REF] Sun | Contour-Shape recognition and retrival based on chain code[END_REF].

The chain code is defined as follows: for one pixel on the boundary of an image or object, it has n neighbors, numbered from 0 to n-1, which are called direction codes. There are 4-direction and 8-direction chain code illustrated in Fig. 4.

/ i i x y A i x x A ∈ = ∑ (3)
( , )

/ i i x y A y y ∈ = i A ∑ (4)
Let r i be the maximum distance from C i to i-direction chain code:

2 ( , ) max ( ( ) (
)

) i i i x y A r x x y ∈ = - + -2 i y (5)
With C i as center and jr i /M as radius, we draw M concentric circles. Let |A ij | be the count of the chain codes with i direction inside jth circle. After normalizing, CCSDE is denoted as

(6) 2 1 ( ) log ( ) M i i ij ij j SE m m m = = -∑
Combined with CCH, the new feature vector is given as 0 0 1 1 ( , ), ,( , ), ,( , )

i i n n h SE h SE h SE - - (7)
If we compare the similarity of two contours c 1 and c 2 , the similarity can defined as

1 2 1 2 1 2 ( ) ( ) 1 ( ) ( ) 1 2 ( ) ( ) 0 min( , ) ( , ) min( , ) max( , ) c c n c c i i i i c c i i i SE SE S c c h h SE SE - = = × ∑ (8)

Color feature extraction

The color features are vital to a scene image but is not sufficient to distinguish similar scenes. To present semantic description better, we combine color feature with the spatial feature. There are a large number of different scene images which have similar structural and spatial features and they are easily classified into the same category using the spatial features only. For instance, Fig. 5 shows two distinctly different scene images of coast and open country but they have similar edge detection images. They both have 3 levels from top to bottom: the former has sky, sea and sands, while the latter has sky, flowers and grassland. The main difference between the two scene images is the color of sands and grassland, instead of the structure and spatial features. Computer system can display over 2 24 colors with 24 bits. We take advantage of this capability to define color data directly as RGB values and eliminate the step of mapping numerical values to locations in a colormap. We define an m-by-n-by-3 matrix to specify truecolor for the scene image. To reduce the number of colors in an image, some functions can be used to approximate the colors in the original image.

Features combination

In our model, the visual words contain two types: RGB color and spatial description. If the image can be easily classified in a category using only color or spatial feature when training the images, the classification is easy done. On the contrary, we take advantage of both color and spatial feature to classify. When using color features, the image is split into block regions. The RGB color features are extracted in every block of each image first and then global RGB color features are formed as visual words. They are vectors of RGB values. After the spatial and color features are extracted, the K-means algorithm is applied to cluster these features. Each cluster center corresponds to a visual word and a visual vocabulary includes K visual words. A set of visual words make up an image, which is similar to that the text is consist of words.

Experiments and results

According to the previous description, in our model, the scene image dataset divides into two parts. One of them is used to extract spatial and color features and train the parameters for the PLSA model, and another is used to test the model we proposed and do experiments to compare with other models. The experiments and the results display as follows.

Datasets and parameters setting

In this paper, we choose Oliva and Torralba (OT) and Fei-Fei and Perona (FP) as the image datasets to conduct our experiment.

The OT dataset has 2688 images and is divided into 8 categories: four nature scenes (328 forest, 374 mountain, 360 coasts and 410 open country) and four artifact ones (308 inside cities, 260 highway, 292 streets and 356 tall buildings).The FP dataset has more than 2688 images but it is only available in gray scale. It includes 13 categories: 8 OT categories plus 174 bedrooms, 151 kitchens, 241 residences, 216 offices and 289 living rooms.

When comparing with [START_REF] Shimazaki | Scene classification using color and structure-based features[END_REF], we use the OT dataset. What's more, we also compare our experiment results with other models based on both datasets.

The values of the latent semantic variables (Z), the number of visual vocabulary (V) in PLSA models and the number of neighbors (K) in KNN classifier are especially important. So, firstly, we choose 100 random scene images from the training set to find the optimal parameters, and the rest of the training images are used to calculate the visual vocabulary and topics of PLSA. However, for the other parameters (the number of block regions etc.), we assign them by referring to the experience or other similar experiments.

Experiments and classification results

As to the classification result, the performance is measured by the average of precision and recall. The precision is defined as the number of images correctly classified in the topic category divides the total number of images in the same category. While the recall is defined as the number of images correctly classified in their categories divides the number of images which are relative to the topic. We use a group of different values of K, V and Z to experiment and analysis the effect and the optimal value. To achieve objective variations, we repeat the experiments 5 times with selecting different random training images and test sets. The mean variation curves are displayed below in Fig. 6.

According to the figures below, we can observe that we get good performance above 75 percent when the number of visual words V=700, the number of latent topics Z=30 and the number of neighbors K=10. We apply our methods to the experiment and use the OT dataset, the same one as [START_REF] Shimazaki | Scene classification using color and structure-based features[END_REF]. The experiment results are comparing in Table 1. The categories of forest, highway, street and tall building have the best result in Ken Shimazaki's experiment. The average is 0.549. However, our proposed model has a better result and almost all 8 categories have a higher mean recall and precision rate. The experiments using color and spatial features respectively and a hybrid of both are done when extracting image features. We use the KNN-SVM classifier we propose. Table 2 displays the performance for both OT and FP dataset by using color and spatial features. After that, the experiments of different classifiers (KNN, SVM and KNN-SVM) are conducted using the method of extracting both color and spatial features we propose when classifying the scene images. Table 3 shows that our classifier is improved the classification performance on OT dataset as well as FP dataset. Other researches on scene classification such as LDA [START_REF] Rasiwasia | Latent Dirichlet Allocation Models for Image Classification[END_REF] and modified PLSA (based on SIFT features) model [START_REF] Bosch | Scene classification via pLSA[END_REF] were proposed. We compare the average performance of our model with them and the results are showed in Table 4.

Conclusion

We have proposed an improved method to classify scene images using PLSA model. When extracting the features of the scene image, the color and spatial features are combined to obtain full data of the image. The spatial features are calculated using the chain code and color features are represented in RGB space. Then we use the PLSA to model and train some images to get the proper parameters of the model. At last, the other images in the dataset are used to test our proposed model. According to our experiments, the results certify that our approach is excellent for scene classification. After the comparison with other models, higher mean precision and recall of our method is achieved, although they are not the best. What's more, the scalability and robustness of our model are also excellent. Our future works are improving the method we proposed to enhance the efficiency and conducting more experiments with different parameters to increase our method's performance.

IFIP IIP2014

China (No. 2012CB326403), the Guangxi Natural Science Foundation (Nos. 2012GXNSFAA053219, 2013GXNSFAA019345) and the "Bagui Scholar" Project Special Funds.

Fig. 1 .

 1 Fig. 1. The coast scene images with different illumination

  And then from P(d,w)=P(d)P(w|d), we get ) is the topic specific distribution and each image is modeled as a mixture of topics P(z|d).The graphical model is represented in Fig.2.

Fig. 2 .Fig. 3 .

 23 Fig. 2. Graphical model of PLSA

Fig. 4 .

 4 Fig. 4. 4-direction and 8-direction chain code

Fig. 5 .

 5 Fig. 5. Similar edge detection images of different scene images

KFig. 6 .

 6 Fig. 6. Performance under variation of different V, Z and K

Table 1 .

 1 Performance comparison between our model and[START_REF] Shimazaki | Scene classification using color and structure-based features[END_REF]

	Categories	Coast	Forest	Mountain	Open country
	Mean Recall and Precision	0.516	0.697	0.354	0.409
	Our model	0.632	0.799	0.532	0.527
	Categories	Highway	Inside city	Street	Tall building
	Mean Recall and Precision	0.670	0.380	0.766	0.656
	Our model	0.783	0.502	0.860	0.837

Table 2 .

 2 Comparison among color, spatial and both features

	# Categories	Color	Spatial	Both (Our model)
	8 OT dataset	0.742	0.683	0.830
	13 FP dataset	0.736	0.622	0.815

Table 3 .

 3 Comparison among KNN, SVM and KNN-SVM classifiers

	# Categories	KNN	SVM	KNN-SVM (Our model)
	8 OT dataset	0.853	0.864	0.871
	13 FP dataset	0.724	0.738	0.743

Table 4 .

 4 Comparison among PLSA (SIFT), LDA and our model

	# training images	128	256	512	1024	1344
	PLSA (SIFT)	0.682	0.753	0.795	0.846	0.859
	LDA	0.772	0.783	0.794	0.810	0.816
	Our model	0.690	0.764	0.798	0.853	0.859

Acknowledge

This work is supported by the National Natural Science Foundation of China (Nos. 61165009, 61262005, 61363035, 61365009), the National Basic Research Program of


Zhuxiao Wang
email: wangzx@ncepu.edu.cn

Jing Guo

Jin Shi
email: shijin@cert.org.cn

Hui He
email: huihe@ncepu.edu.cn

Ying Zhang
email: yingzhang@ncepu.edu.cn

Hui Peng
email: penghui@bisu.edu.cn

Guanhua Tian
email: guanhua.tian@ia.ac.cn

A Formal Model for Attack Mutation Using Dynamic Description Logics

Keywords: Intrusion Detection/Prevention Systems, Multi-protocol Evasions, Advanced Evasion Techniques, Knowledge Representation and Reasoning, Dynamic Description Logics

All currently available Network-based Intrusion Detection Systems (NIDS) rely upon passive protocol analysis which is fundamentally flawed as an attack can evade detection by exploiting ambiguities in the traffic stream as seen by the NIDS. We observe that different attack variations can be derived from the original attack using simple transformations. This paper proposes a semantic model for attack mutation based on dynamic description logics (DDL(X)), extensions of description logics (DLs) with a dynamic dimension, and explores the possibility of using DDL(X) as a basis for evasion composition. The attack mutation model describes all the possible transformations and how they can be applied to the original attack to generate a large number of attack variations. Furthermore, this paper presents a heuristics planning algorithm for the automation of evasion composition at the functional level based on DDL(X). Our approach employs classical DL-TBoxes to capture the constraints of the domain, DL-ABoxes to present the attack, and DLformulas to encode the objective sequence of packets respectively. In such a way, the evasion composition problem is solved by a decidable tableau procedure. The preliminary results certify the potential of the approach.

INTRODUCTION

A weakness of most currently available Network-based Intrusion Detection Systems (NIDS) that rely upon passive protocol analysis is their inability to recognize an attack that evades detection by exploiting ambiguities in the traffic stream as seen by the NIDS. Exploitable ambiguities may arise in different ways: [START_REF] Ptacek | Insertion, evasion, and denial of service: Eluding network intrusion detection[END_REF] The NIDS may lack complete analysis for the full range of behavior allowed by a particular protocol.

(2) Without detailed knowledge of the victim end-system's protocol implementation, the NIDS may be unable to determine how the victim will treat a given sequence of packets if different implementations interpret the same stream of packets in different ways. [START_REF] Niemi | Protect Against Advanced Evasion Techniques[END_REF] Without detailed knowledge of the network topology between the NIDS and the victim end-system, the NIDS may be unable to determine whether a given packet will even be seen by the end-system [START_REF] Ptacek | Insertion, evasion, and denial of service: Eluding network intrusion detection[END_REF][START_REF] Handley | Network Intrusion Detection: Evasion, Traffic Normalization, and End-to-End Protocol Semantics[END_REF]. Advanced Evasion Techniques (AET's) is a lately established term in network security industry referring to a set of non-trivial and expensive means to bypass NIDS in order to deliver an exploit, attack, or other form of malware to a target network or system, without detection. Advanced evasion techniques can be identified according to certain underlying principles: (1) Delivered in a highly liberal way; (2) Employ rarely used protocol properties; (3) Use of unusual combinations; (4) Craft network traffic that disregards strict protocol specifcations; [START_REF] Wang | Failure Diagnosis of Internetware Systems Using Dynamic Description Logic[END_REF] Exploit the technical and inspection limitations of security devices: memory capacity, performance optimization, design flaws [START_REF] Niemi | Protect Against Advanced Evasion Techniques[END_REF].

Since we are interested in testing the ability of a NIDS to properly identify real intrusions, we need a way to ensure that executing each of mutants generated by advanced evasion techniques against the vulnerable application, we are going to obtain the same effect as executing the original exploit script. In this article, we define an attack mutation model that describes all the possible transformations and how they can be applied to the original attack to generate a large number of attack variations. The semantic model for attack mutation is based on dynamic description logics (DDL(X)), extensions of description logics (DLs) with a dynamic dimension. The attack mutation model is self explanatory. Given an original attack, the attack mutation model can provide a proof that a sequence of transformations used for obfuscation is a real attack. Developers can use the model to analyze attacks and to identify the exact transformation that their NIDS fails to handle. The attack mutation model is exhaustive, capable of generating all attack variations from a known base attack using a set of rules. All the mutation techniques are individually sound, and also that any possible composition of them is sound. So the model is sound, generating only instances that implement the original attack.

In this article, we present an approach for automatic evasion method plans based on dynamic description logics [START_REF] Chang | A Family of Dynamic Description Logics for Representing and Reasoning About Action[END_REF][START_REF] Wang | Failure Diagnosis of Internetware Systems Using Dynamic Description Logic[END_REF][START_REF] Wang | An architecture description language based on dynamic description logics[END_REF] named DDL(X), extensions of DLs [START_REF] Baader | The description logic handbook: theory, implementation, and applications[END_REF] with a dynamic dimension [START_REF] Artale | A temporal description logic for reasoning about actions and plans[END_REF][START_REF] Baader | Integrating description logics and action formalisms: First results[END_REF]. Our approach used classical DL-TBoxes to capture the domain constraints, DL-formulas to encode the objective sequence of packets, and DL-ABoxes to be a special representation of the attack that provides to the underlying mutation engine the mechanism to manipulate the attack, respectively. Actions in DDL(X) were used to abstract the functionalities of the available evasion methods. In such a way, the automatic evasion plans can be reduced to formula satisfiability checking in DDL(X) and solved by a decidable tableau procedure.

In the following sections, we firstly illustrate how variants of a real exploit can be derived from the original exploit script in Section 2. In Section 3, we demonstrate the descriptions of evasion techniques can be formalized as actions in DDL(X), and formalize the notion of the attack mutation model using dynamic description logics. Afterwards, in Section 4 we detail the problem of evasion composition can be solved by reasoning about actions in DDL(X) and present a heuristic planning algorithm for automated composition of evasions. Finally, we summarize the paper in Section 5.

APPROACH OVERVIEW

First of all, we need an instance of the exploit script that we want to mutate. The base instance is then used to derive another attack instance by repeatedly applying single step transformations. Our example vulnerability is a published buffer overflow in a commonly used Windows XP SP2 host (MSRPC Server Service Vulnerability, CVE-2008-4250 in www.cve.mitre.org.); exploiting the overflow may allow arbitrary code execution. The exploit causes the overflow by providing a crafted RPC request that triggers the overflow during path canonicalization. We call this exploit E MSRPC .

Given an exploit E, Trace(E) denotes a sequence of packets from a NIDS perspective and the function Post(E) to be the post-conditions of the execution of E against a target system. Consider the following twelve transformation rules in Table 1 that can be applied to an existing operational description of how a vulnerability is exploited to generate a new different version of the same exploit.

We call these rules semantics preserving because they do not alter the semantics of E, i.e. the transformation does not affect the results of the execution of the exploit. If E is an instance of the E MSRPC attack, then by using rules in Tab.1 it is possible to derive the conclusion that the E' is also an instance of E MSRPC and the instance E' contains the necessary data for a successful E MSRPC attack. For example, we apply r 1 on E to send extra NetBIOS packets to break the packet flow. Then, we apply r 9 and change the order of TCP segments. In dynamic description logic terminology, starting with an exploit E, we can successively apply a set of transformations T ={r 0 , r 1 , …, r n } to create a complex mutant exploit E'.

To formalize the notion of the semantics preserving transformation, we get:

Trace(E) ≠ Trace(E') Post(E) = Post(E')
The first condition requires that the transformation manifests itself as a change in the sequence of packets. The second condition requires that the transformation preserves the attack post-conditions. While the E and the E' might look different from a NIDS perspective, from the attacker point of view the E' generated by evasion technique is still an "effective" attack (i.e., that executing the E' against the vulnerable application, we are going to obtain the same effect as executing the original exploit script E). Intuitively speaking, one can infer E' from E, and vice versa. Change the order of IPv4 fragments

In our description frameworks for evasions, functional descriptions are essentially the state-based and use at least pre-state and post-state constraints to characterize intended executions of an evasion. Many evasion techniques can be combined. When evasion techniques are combined at different levels, NIDSs that detect each separate evasion technique often fail at detecting some permutations. First, we apply the transformation rules in a Breadth-First order: we first apply application level rules because they are independent against TCP-level and IP-level rules, then we segment each instance into small pieces, change the order of TCP segments we get, and send chaff IPv4 packets interleaved with normal packets(Figure 1). Second, we prune away some of the derivation branches to decrease the number of instances. In some cases, only a subset of possible factors evades NIDS detection; adding additional evasion measures to an evasive attack may cause the NIDS to detect it.

A FORMAL MODEL FOR ATTACK MUTATION

A DDL(X) model is a tuple M = (W, T, Δ, I), where, W is a set of states;

T : N A →2 W×W is a function mapping action names into binary relations on W;

Δ is a non-empty domain; Definition 1 (Atomic evasions) An atomic evasion is a tuple t=<Pre, Effects>, where Pre is a finite set of formulas in DDL(ALCO) specifying the preconditions for the execution of t; and Effects is a finite set of assertions or their negation in ALCO, which is the facts holding in the newly-reached world by the evasion's execution.

Composite evasions are constructed from atomic evasions with the help of classic constructors in dynamic description logics [START_REF] Chang | A Family of Dynamic Description Logics for Representing and Reasoning About Action[END_REF].

Below we formally define a mutation model for an exploit as well as some reasoning tasks and the planning problem.

Definition 2 (Mutation Model of an Exploit)

Let E be an instance of an exploit and T be a set of sound inference rules with respect to E.

A mutation model of E is a DDL(X) model (W, T, Δ, I).

Such a model enables derivation of new exploits by applying the inference rules (like evasion methods) on already known exploits. For an attack E, we envision the attack mutation model that, with respect to a set of rules, is sound: derives only sequences of packets that implement E; complete: can derive any sequence of packets that implements E; and decidable: given a sequence of packets, there is an algorithm that determines whether or not the sequence is derived from the already known exploit.

Definition 3 (NIDS View)

Let N be a NIDS. N's view with respect to an exploit E, denoted S E N , is the set of sequences of packets that N recognizes as E. Given a NIDS and an instance of an exploit E, the basic reasoning task for DDL(X) is to find an instance of E that evades the NIDS. Definition 4 (Reasoning Task 1) Let (W, T, Δ, I) be an attack mutation model of E, and N be a NIDS. Let S E N be the view of N with respect to E. The reasoning problem is to find a sequence of packets S that is derivable from E, but is not in S E N . More formally, find S  S E N such that M ⊨ S and M ⊨ E. Given an instance of an exploit E and a sequence of packets S, another important reasoning task for DDL(X) is to determine whether S is an instance of E. The last inference problem we will investigate is the planning problem. Given a goal statement (i.e. a sequence of packets) and a set of actions(i.e. evasion methods), the planning problem is to find an action sequence that will lead from the initial state(i.e. an already known exploit) to states in which the goal will hold. With DDL(X), we define the plans as follows.

Definition 6 (The Planning Problem) An action sequence r 1 , ..., r n is a plan for a goal S w.r.t. M = (W, T, Δ, I) if and only if (i) the sequence-action r 1 ;...; r n is executable on states described by E and (ii) S is a consequence of applying r 1 ;...; r n on states described by E.

It is intuitive to model evasions by actions in DDL(ALCO). As demonstrated in this section, the functionalities of evasions can be semantically transformed into actions in DDL(ALCO) by a proper domain ontology (TBox). All kinds of reasoning tasks concerning the functionalities of evasions thus can be reduced to the reasoning about actions in DDL(ALCO), which are the topic of the next section.

AN EFFICIENT ALGORITHM FOR EVASION COMPOSITION

In this section we demonstrate that the problem of evasion composition can be reduced to satisfiability checking of formulas in DDL(ALCO) and then be solved by a decidable procedure after the transformation process that transforms evasions to actions in DDL(ALCO). Afterwards, we propose a heuristic planning algorithm for automated composition of evasions. The algorithm achieves good balance between computational performance and accuracy.

When facing a problem of evasion composition, we firstly collect the relevant evasions and transform these preexisting evasions to atomic actions by constructing the specification of the domain.

Let us analyze the following two formulas:

( [ (α 1 …α n )* ] Π  Con j(E) )  < plan > true,
where Conj(E) represents the conjunction of all the elements of the set E, Π the formula  1 n i (Con j(P i )  <α i > true) and P i the precondition of action r i for each i: 1≤i≤n.

The whole formula above states that the action sequence plan is executable on states described by E. The above formula is valid if its negation (labeled as Eq. 1) is not satifiable:

[ (α 1 …α n )* ] Π  Con j(E)  < plan > true (1) Con j(E)  [plan]S .
This formula indicates that the goal S is a consequence of applying the action sequence plan on world states described by E. Similarly, the above formula is valid if its negation (labeled as Eq. 2) is not satifiable:

Con j(E)   [plan]S (2)
Given a goal statement S and a set of actions ∑ = {r 1 , r 2 ,…, r m }, Algorithm 1 shows how to produce a plan (i.e. an action sequence) that will lead from the initial state to states in which the goal S will hold. The heuristics algorithm ActionPlan() travels the possible world states in a Breadth-First manner and terminates at a successful plan, or failure after the algorithm attempts nearly exhaustively. However, instead of travelling all possible world states, the heuristic selects world states more likely to produce an successful plan than other world states. It is selective at each decision point (Line.11), picking world states that are more likely to produce solutions.

Algorithm 1 ActionPlan(E, T, ∑, S, plan)

Input: initial ABox E; TBox T; the set ∑= {r 1 , r 2 ,…, r m } of available actions; and formula S, seen as a goal statement.

Output

SUMMARY

In this article, we aim at providing an effective framework for the composition of evasion techniques to test the quality of intrusion detection signatures. Our approach supports multiple evasion techniques and allows the developer of the test to compose these techniques to achieve a wide range of attack mutations. We define an attack mutation model that describes all the possible transformations and how they can be applied to the original attack to generate a large number of attack variations. The attack mutation model is based on DDL(X), extensions of description logics (DLs) with a dynamic dimension. In particular, we proposed a heuristic planning algorithm for automated composition of evasion methods by a reduction to the formula satisfiability checking in DDL(X) and a selection of world states more likely to produce a successful plan. The functionalities of the evasion methods are abstracted by actions in DDL(X), while the domain constraints, exploits, and the objective sequence of packets are encoded in TBoxes, ABoxes and DL-formulas, respectively.

Then the problem of evasion composition can be reduced to formula satisfiability checking in DDL(X) and solved by a decision procedure. Afterwards, instead of travelling all possible world states in a Breadth-First manner, a heuristic planning algorithm selects world states more likely to produce an successful plan than other world states. Our approach has several important advantages: the attack mutation model permits the application of analytical methods for deriving sound evasion composition; executing each of mutants generated by our approach against the vulnerable application, we are going to obtain the same effect as executing the original exploit script.

Fig. 1 .

 1 Fig. 1. An Illustration of a Sample Multi-protocol Evasion I is a function which associates with each state w  W a description logic interpretation I(w) =< Δ, • I(w) >, where the mapping • I(w) assigns each concept to a subset of Δ, each role to a subset of Δ×Δ, and each individual to an element of Δ.

Definition 5 (

 5 Reasoning Task 2) Let M = (W, T, Δ, I) be an attack mutation model of E and S be a sequence of packets. S is an instance of E if and only if there exists a model M = (W, T, Δ, I) and a state w  W such that M ⊨E and (M, w)⊨S.

Table 1 .

 1 Atomic evasions

		Name	Description
	NetBIOS	NetBIOS chaff(r 1 )	Send extra NetBIOS packets to break the packet flow
		NetBIOS initial chaff(r 2 )	Send chaff NetBIOS packets when establishing the NetBIOS connection
	SMB	SMB filename obfuscation(r 3 )	Obfuscate the tree name used in the SMB NT Create AndX method
		SMB WriteAndX padding(r 4 )	Insert extra padding between the WriteAndX header and payload
	MSRPC	MSRPC request segmentation(r 5 )	Set the maximum number of bytes written in a single MSRPC fragment
		MSRPC NDR modifications(r 6 )	Set NDR types not related to endianness
	TCP	TCP Chaff(r 7 )	Send chaff TCP segments to baffle inspection
		TCP timestamp option settings(r 8 )	Set initial TCP timestamp option settings
		TCP segment order(r 9 )	Change the order of TCP segments
	IP	IPv4 chaff(r 10 )	Send chaff IPv4 packets interleaved with normal packets
		IPv4 fragmentation(r 11 )	Fragment IPv4 packets to given size
		IPv4 fragment order(r 12 )

Algorithm 2 EvaluatePlan(E, T, goal, srcplan, dstplan) Input: initial

  ABox E; TBox T; formula goal, seen as a goal statement; srcplan, a plan for left-hand side of comparison; dstplan, a plan for right-hand side of comparison.Output: an integer greater than, equal to, or less than 0, if the number of primitive formula achieved through dstplan is greater than, equal to, or less than the number of primitive formula achieved through srcplan, respectively.

	16. end if
	17. end while
	18. return nil.
	End
	Begin
	1. srccount←0;
	2. dstcount←0;
	3.	rewrite formula goal to a disjunctive normal form formula
	G= 0 d i ( 0 m j ϕ ij );
	4. for each primitive formula ϕ ij do
	5. if Con j(E)  [srcplan] ϕ ij is unsatifiable then
	6.	srccount←srccount+1;
	7. end if
	8. if Con j(E)  [dstplan] ϕ ij is unsatifiable then
	9.	dstcount←dstcount+1;
	10. end if
	11. end for
	12. return dstcount-srccount;
	End	: an successful plan or nil as failure;
	Begin
	1. initialize plan with an empty list;
	2. Initialize queue QueOfPlans with plan;
	3. while (QueOfPlans is unempty) do
	4. remove head of QueOfPlans and set it to plan;
	5. if (Eq. 2 is unsatifiable) then
	6. return plan as a successful plan;
	7. else
	8. for each r i ∑ do
	9. 10.	newplan ← < plan, r i >;// appending r i to the rear of the list plan. if([ (r 1 …r k )* ] ΠCon j(E)  < newplan > true is unsatifiable)
	then
	11.	if(EvaluatePlan(E, T, S, plan, newplan) ≥ 0) then
	12.	queue QueOfPlans with newplan;
	13.	end if
	14.	end if
	15. end for

ACKNOWLEDGEMENTS

This work is supported by the National Natural Science Foundation of China (No. 61300132, No.61363058, No.61372182, No.61305056, No.61303172), the Specialized Research Fund for the Doctoral Program of Higher Education (No.20120036120003), the Fundamental Research Funds for the Central Universities (No.12ZP09, No.13ZP010) and the Scientific Plan Project of Beijing Municipal Commission of Education(No.SQKM201420131001).