Abstract

An artificial bioindicator system is developed in order to solve a network intrusion detection problem. The system, inspired by an ecological approach to biological immune systems, evolves a population of agents that learn to survive in their environment. An adaptation process allows the transformation of the agent population into a bioindicator that is capable of reacting to system anomalies. Two characteristics stand out in our proposal. On the one hand, it is able to discover new, previously unseen attacks, and on the other hand, contrary to most of the existing systems for network intrusion detection, it does not need any previous training. We experimentally compare our proposal with three state-of-the-art algorithms and show that it outperforms the competing approaches on widely used benchmark data.

1 Introduction

This paper proposes a novel approach for tackling the network intrusion detection problem. The presented model is based on the effects of environmental changes or perturbations on a population of highly sensitive evolutionary agents. The agent population resides in a system (the environment) that is being monitored, and the agent population serves as a bioindicator [30] for detecting abnormal system behavior. In other words, the quantification of the perturbation effect on the population of agents is used for detecting abnormalities in the monitored system. For the purpose of this work, the monitored system is a computer network.

One possible solution for the problem of anomaly detection consists in classifying the behaviors of the monitored system. In some cases, the classes are predefined and do not change over time. In more complex cases, classes may not be defined a priori, and may change over time. The detection of intrusion in data networks, which is the aim of network intrusion detection systems (NIDSs), is one of the most complex scenarios for anomaly detection. In this domain, the classifying algorithm must deal with at least two fundamental classes: normal traffic and intrusive traffic. These classes are not static, as they change due to the usual variation in the behavior of system users or the presence of new or unknown attacks. Hence, this scenario has been chosen to test the capacities of the classifier proposed in this paper.

1.1 Relation to Existing Works

The network intrusion detection problem has been tackled by a wide range of approaches. The interested reader may find a complete taxonomy in [41]. In the following we will give a brief survey.

The most common approach for the network intrusion detection problem in industry is signature detection [7]. This kind of technique uses a pre-coded representation of a set of intrusions. Moreover, the system is characterized by its lack of learning capacities. Signature detection is based on a static set of intrusion detection rules, which are generated beforehand by the use of expert knowledge. The advantage of this approach is its ability to tag attacks, in other words, the skill of the system for identifying attack classes. Its disadvantages are its dependence on signature updates [48] and its inability to detect new attacks [22]. Signature detection is used, for example, by Cisco, SourceFire, Checkpoint, and Symantec as well as in other commercial solutions.

Concerning machine learning techniques, a variety of different methods, including neural networks [52], support vector machines [40], genetic programming [44], fuzzy logic [57], Bayesian networks [3], and hidden Markov models [9], among others, have been used to address the network intrusion detection problem. A good overview of recent developments may be obtained from [56, 62, 35, 7, 20, 41]. The disadvantage of these kinds of techniques is that they require separate processes of training under controlled conditions, and high-quality labeled training data sets. Moreover, they are unable to detect unknown attacks or adapt to scenario changes, because of their static initial training.

Finally, the intrusion detection problem has also been tackled by techniques based on statistics. These approaches make use of the network traffic activity for profiling the system [63]. Various clustering techniques have also been applied to this problem, such as, for example, self-organizing maps [50].

1.2 Relation to Artificial Immune Systems

The system proposed in this work is closely related to artificial immune systems (AISs). Like that of their natural counterparts, the task of artificial immune systems is the detection and elimination of components that are not own, harmful, or destabilizing to the system in which they reside. AIS research began with Farmer, Packard, and Perelson's study [17] in 1986. One of the first—and most explored—theories is negative selection (NS) [19], which is based on the maturation of T lymphocytes to produce immunological tolerance (see [25]) and on models based on the Jerne immune network [14]. Among other applications, techniques based on negative selection have been successfully used for intrusion detection problems [33, 46]. For a recent survey of negative selection techniques, the interested reader is referred to [31]. Other examples of approaches belonging to the paradigm of AIS include the clonal selection algorithm (CSA) [12], which is based on a proliferation of detectors capable of detecting antigens and exploring them in order to enhance affinity by means of somatic hypermutation. Examples of systems based on the CSA include [10, 23]. Moreover, an application to intrusion detection was proposed in [34]. Other aspects of natural immune systems are used in the immune network algorithm (INA) [13], such as, for example, the aiNet [14] and danger theory algorithms [45]. The inspiration for these systems is the immune activation process. They were applied to intrusion detection in [21, 2]. In general, all these approaches aim at the identification and the elimination of threats, while maintaining failure tolerance and adaptability.

However, the field of AISs in engineering is still young, and several debates are still ongoing. Bersini [5], for example, argues that approaches based on the traditional conception of the immune system as a purely defensive entity, which is adopted in the approaches outlined before, may be incorrect. In any case, the importance of artificial immune systems for engineering may be found in the concept of endogenous double plasticity [11], which states that a system adjusts structurally during its functioning in a continuous and plastic manner, integrating new elements and discarding old ones, with the changes being controlled by its internal dynamics. Such a behavior may be achieved by simple heuristics such as compensating for weak elements, maintaining diversity, and eliminating redundancy, that is, maintaining a balance through ecological mechanisms. This is exactly the aim of the artificial bioindicator system presented in this paper.

1.3 Contribution of this Work

The artificial bioindicator system that is the subject of this work presents an ecological approach to network intrusion detection, assuming that an agent population will evolve learning skills in order to survive in the system in which it resides by adapting to the characteristics of the system. This can be seen as a structural modification of the agent population. This adaptation process allows the transformation of the agent population into a bioindicator capable of reacting to a system anomaly. This ecological approach is present, for example, in Varela's constructivist vision of the biological immune system [59, 6, 58], which emphasizes self-affirmation and homeostatic potential.

Our artificial bioindicator system aims to achieve substantial operational advantages over existing techniques for network intrusion detection, which—as outlined before—are mostly limited to classifiers relying on the use of recognized attack signatures or patterns [47]. The main drawback of existing approaches is that they require continual supervised updating processes in order to remain useful. Otherwise, they rapidly become ineffective against unknown attacks [22]. Although there are existing approaches that provide a solution to the problem of detecting novel attacks, they still come with the drawback of yielding a significant number of false positives [16]. Our aim with the artificial bioindicator system is to solve these problems by taking advantage of the anomaly detection abilities of artificial immune systems. Finally, we would like to point out that our proposal may also be classified as a so-called concept drift approach [38], because artificial bioindicator systems constitute an active learning procedure, learning and forgetting continuously.

The remainder of this paper is organized as follows. We describe the artificial bioindicator system in Section 2, and the classifier calibration in Section 3. Section 4 introduces the testing procedure. Finally, the results are presented in Section 5, and conclusions in Section 6.

2 The Artificial Bioindicator System

In the following section we describe the artificial bioindicator system that we propose for network intrusion detection. Henceforth, we refer to the system as a flow anomaly detection system (FADS). The proposed FADS analyzes the behavior of the network by means of a population of agents through continuous observation of a predefined set of external parameters. These parameters characterize the activity of the monitored system. The system is inspired by Varela's view of the functioning of the immune system [59, 6, 58]. In that view, the immune system is a symbiont that protects its host, maintaining its equilibrium through ecological mechanisms. In the context of this paper, the host corresponds to the monitored system, while the population of agents is the symbiont that continuously adapts to the host. This is called a structural coupling.

The FADS that was developed in this work has been implemented in NetLogo [51]. The system can technically be described as follows. It consists of a bidimensional w × h grid topology T that is inhabited by a set A of evolutionary agents. T is also called the environment. Positions in T are denoted by (i, j) where 1 ≤ iw and 1 ≤ jh. Each agent aA is characterized as follows:

  • 1. 

    An agent a has a fixed position (ia, ja) in T.

  • 2. 

    Moreover, a has a variable energy value ea ≥ 0.

  • 3. 

    Finally, a has a genetic code πa, which is a permutation of the integer values in {0,…, n}. Each integer refers to one of n + 1 particle types, which may be found by moving through the environment T.

Apart from the set A of agents, the environment is also inhabited by a set P of moving particles. Each particle pP is characterized as follows:
  • 1. 

    A particle pP has a current position (ip, jp) in T.

  • 2. 

    Moreover, p has a variable energy value ep ≥ 0.

  • 3. 

    Finally, p also has a particle type tp ∈ {0,…, n}.

Note that particle types (which are represented as integer values) correspond to features of the monitored system. In the context of network intrusion detection, particle types correspond to data packet types. The creation of a particle p of type tp in the system is determined by the presence of the corresponding feature in the current activity of the monitored system. Thus, the FADS receives as inputs binary vectors representing observations of the monitored system. As previously mentioned, agents have fixed positions in T. In contrast, particles move at each iteration from their current position to one position to the left. In those cases in which particles find themselves in a position of an agent, they interact in the way described below.

The pseudocode for the simulation of the system is presented in the following algorithm, and a schematic representation in Figure 1 (the procedure for finding well-working parameter values is described in detail in Sections 3 and 4):

graphic

Figure 1. 

Graphical presentation of the topology used. The particles pP (the black squares) enter the topology T on the right side, and leave it on the left. The agents aA (gray dots) have fixed positions assigned to them at birth.

Figure 1. 

Graphical presentation of the topology used. The particles pP (the black squares) enter the topology T on the right side, and leave it on the left. The agents aA (gray dots) have fixed positions assigned to them at birth.

The input of the system consists of a set A of agents with random positions in T and a random genetic code. The initial energy value of each agent is set to a value einit > 0, which is a parameter of the system. Moreover, A(t) refers to the agent population at iteration t, and P(t) refers to the particles in the system at iteration t. At each iteration, the following steps are performed. First, in lines 5–8 all existing particles move one position to the left in T (see (ii) in Figure 1). When a particle leaves the environment, it is removed from the set P of currently existing particles (see line 7, and (iv) in Figure 1). Secondly, new particles are generated in the function GenerateNewParticles() (line 9). This is done by a continuous observation of the system that is to be monitored. Each of the n + 1 observable characteristics, in the case of detection, is translated into a new particle p of type tp corresponding to the observed characteristic. That is, if the observed characteristic is c ∈ {0,…, n}, then a new particle p of type tp = c is generated. The position of this new particle is chosen randomly from {(w, 1),…, (w, h)}, that is, a random position from the rightmost column of T. This is done under the restriction that no two particles may reside in the same position of T (see (i) in Figure 1).

Finally, agent actions are performed in lines 10–19 of Algorithm 1. First, the metabolism of each agent causes an energy consumption of m units at each iteration (see line 11). Moreover, if an agent a finds itself in the same position as a particle p (line 12, and (iii) in Figure 1), both the energy value of the agent and the energy value of the particle are affected. The change in the energy value of the agent depends on the type tp of particle p. More specifically, the function AdaptEnergyOfAgent() produces a change in agent a's energy value as follows. First, the nutritional value ζ of particle p for agent a is determined as follows:
formula
where ϕ is the maximal nutritional value that a particle may provide, and ϵ is a parameter that determines the linear nutritional loss. Both ϕ and ϵ are parameters of the system. To conclude, k corresponds to the position of tp (the type of particle p) in the genetic code πa of agent a. Note that the smaller k, the higher is the nutritional value provided by particle p to agent a. Finally, the nutritional value ζ is added to the agent's energy value, ea. In the function AdaptEnergyOfParticle() the energy value of particle p is reduced by exactly one unit (see line 15). When the energy value of an agent or a particle has reached (or fallen below) zero, the agent or particle is considered to be exhausted and is removed from the topology (see lines 14 and 16). The initial energy value of each particle is set to a value vinit > 0, which is a parameter of the system.

The last action of each iteration concerns the creation (birth) of new agents. An agent a engages in asexual reproduction when its energy value ea is greater than or equal to some threshold value r > 0. If this is the case, a new agent a′ is generated and added to A. The position of the new agent is chosen randomly in the von Neumann neighborhood of agent a in T.1 Furthermore, the genetic code πa of the new agent a′ results from a copy of the genetic code πa of agent a in which two random positions are swapped. Finally, the energy value of the agent a′ is set to an initial value ebirth, which is a parameter of the algorithm. Note that, for the purpose of energy preservation in the system level, the energy value ea of agent a is reduced by the sum of ebirth and eproc, as that is the total energy consumed in the process.

The procedure DetectAnomaly(timesOutOfBound,t) implements the classification stage of the algorithm. The pseudocode can be found in the following algorithm:

graphic

This function uses the real-time status of the agent population in order to decide whether to launch an attack alert. First, in line 2 the average number of agents per iteration over the last ψ iterations is computed. This number is denoted by SA. Then, in line 3 the average energy of the agent population over the last ρ iterations is calculated. The notation for this number is Se. In this context, remember that the current iteration is labeled t. In line 4, the algorithm checks if either |A(t)| or e(t) is out of bounds. The definition of being out of bounds is as follows. If |A(t)| surpasses (or falls below) SA by more than dA, where dA is a fixed threshold value, |A(t)| is said to be out of bounds. Similarly, whens e(t) surpasses (or falls below) Se by more than de, where de is also a fixed threshold value, e(t) is said to be out of bounds. Note that being out of bounds indicates a perturbation of the normal system behavior. However, during experimentation we noticed that launching an attack alert whenever the system behavior is out of bounds would generate too many false positives. Therefore, it was decided to introduce the variable timesOutOfBound, which counts for how many consecutive iterations the system behavior has been out of bounds. If the value of this variable surpasses a threshold u, the system launches an alert by calling the function LaunchAlert() (see line 10). However, this is only done if the current iteration t is greater than tstart. The purpose of this is to give the agent population tstart iterations at the start of the simulation in order to stabilize. In other words, the variable timesOutOfBound implements an alert delay with the intention of avoiding unnecessary false positives. After initial experiments, a value of tstart = 2000 was chosen throughout the paper. The graphic of Figure 2 illustrates the working of procedure DetectAnomaly(timesOutOfBound,t).

Figure 2. 

Illustration of the behavior of the agent population. The graphic shows the evolution of the average number of agents per iteration (labeled |A(t)|), the evolution of the average energy of the agent population e(t), and the evolution of the threshold lines (SAdA), (SA + dA), (Sede), and (Se + de) as determined by the threshold values dA and de. Notice that the system first needs to stabilize. Then, after about two-thirds of the timeline, an attack occurs. This causes |A(t)| and e(t) to drop below the threshold lines.

Figure 2. 

Illustration of the behavior of the agent population. The graphic shows the evolution of the average number of agents per iteration (labeled |A(t)|), the evolution of the average energy of the agent population e(t), and the evolution of the threshold lines (SAdA), (SA + dA), (Sede), and (Se + de) as determined by the threshold values dA and de. Notice that the system first needs to stabilize. Then, after about two-thirds of the timeline, an attack occurs. This causes |A(t)| and e(t) to drop below the threshold lines.

At this point, it is important to note that particle types do not correspond to possible classes of the classification problem of differentiating between normal and pathological network traffic. Particle types rather correspond to behavioral characteristics of the host system under observation. The agent population may, or may not, be able to adapt itself to certain particle types. Particles to which the agents are adapted are beneficial and increase the energy of the agents. Particles to which the agents are not adapted will cause a decline in energy or even death. This behavior allows the agent population to observe the presence of an anomaly. In this way, the population of agents can be seen as a biomarker with the ability to detect an organic disorder (corresponding to Varela's view of the immune system). With respect to the metaphor of immune systems, agents might be seen as lymphocytes, but not in the traditional way in which these agents are seeking particular patterns of hazards (pathogens) to mark or to destroy. This paper understands lymphocytes as ecological agents that make adjustments to suit their environment (host), aiming at reestablishing the balance of the system and thereby allowing the survival of the entire system (host and symbiont).

3 System Calibration

A classifier based on the FADS presented in the previous section requires a set of sensitive agents, that is, a set of agents that rapidly change their state in a quantitative manner when the incoming particle flow changes. Remember that such a change of the state of the agent population can be detected via an observation of the evolution of the variables |A(t)| (the current number of agents) and e(t) (the current average energy level of the agents). It is reasonable to assume that the population of agents attains sensitivity to changes of the incoming particle flow (representing environmental changes) from an adjustment of the system parameters. Moreover, a parameter setting is favorable for our purposes if it keeps the agents in a state of stress and if it is highly unfavorable for the consolidation of the agent population.

For the purpose of measuring the quality of different parameter settings, the so-called sensitivity indicator is defined in the following as a measure of the reaction of |A(t)| and e(t) to variations in the incoming flow of particles. More specifically, this measure is defined as the ratio between the number of detected flow variations (Δd) and the real number of flow variations (#Δ):
formula
Here, Δd is the number of LaunchAlert() events at times when a flow variation really occurred (true positives). Moreover, note that in the context of the experiments concerning the system calibration, the parameter u was set to zero in order to achieve increased system sensitivity.2
As a second measure for the quality of a parameter setting we will use the population diversity ψ(t) of the agents at iteration t, which is defined as follows:
formula

In the following, we first summarize the set of system parameters that were subject to system calibration. Afterwards, we introduce the function (henceforth called the flow generation function) that was used to generate the incoming flow of particles, as well as changes in that flow.

3.1 Summary of System Parameters

The developed FADS constitutes a system comprising an agent population that adapts to the normality of a monitored system, that is, the incoming flow of particles. Due to the interest in achieving an agent population reacting in a sensitive and quantifiable manner to disorders caused by abnormalities (changes in the particle flow), the following eight system parameters are subject to system calibration:

  • • 

    ebirth: This is the initial energy level assigned to a newborn agent (created by the function CreateNewAgent() of Algorithm 1). Due to energy conservation laws (which are usual in artificial life systems [53]), this energy is discounted from the progenitor.

  • • 

    einit: This is the energy level of the initial set of agents. Note that we differentiate between ebirth and einit, because agents from the initial agent population need a higher energy level due to their random genetic code, which makes it more difficult for them to survive.

  • • 

    eproc: The value of this parameter corresponds to the energy consumed by an agent for producing an offspring. Again, due to energy conservation laws, this quantity is taken away from an agent when engaging in reproduction. Therefore, an agent loses the sum of ebirth and eproc for producing an offspring.

  • • 

    r: The value of this parameter determines the minimum energy needed by an agent to reproduce.

  • • 

    vinit: The initial energy level of a particle. Note that its value also determines the maximal number of times a particle may interact with an agent.

  • • 

    m: This parameter determines the baseline energy consumption of each agent per iteration of the system (see line 11 of Algorithm 1).

  • • 

    ϕ: This parameter corresponds to the maximal nutritional value that is needed to calculate the nutritional value ζ of a particle for the agent with which it interacts. In other words, its value contributes to determining the energy delivered by a particle to an agent (see also Equation 1).

  • • 

    ϵ: This is another parameter that is needed for the calculation of the nutritional value ζ of a particle for the agent with which it interacts. This parameter determines the energy loss of agents in the event of a lack of genetic suitability with respect to the particles with which they interact.

3.2 Flow Generation Function

For the purpose of calibrating the system parameters, we developed the following flow generation function based on a normal distribution. We did not use real network traffic, in order to avoid the overfitting of the system to a specific application. The flow generation function can be used for producing 26 different particle flows (as described below), thus allowing for the selection of system parameters in order to achieve a sensitive system.

As previously outlined, the input to the system are particles of n + 1 different types. The flow generation function is a normal distribution defined over a binary vector of length n + 1. The positions in this vector are (at each iteration) set to zero. Then, the flow generation function is sampled until the binary vector contains exactly #p positions that are set to one. A value for #p, which is the number of input particles per iteration, must be determined before the start of the simulation. Note that, when a position i in the binary vector is set to one, this means that a particle of type i is generated and fed into the system. We used a setting of n = 49 for the calibration of the system, which corresponds to a number of 50 different packet types. This setting was inspired by works from the literature on this topic, which used around 50 different packet types. At each iteration, a binary vector of length n + 1 is generated, based on the parameters described in the list below. For each parameter we provide the base setting and the set of actions that may be performed for the definition of a different normal distribution, resulting in a variation of the input of particles into the system (henceforth called a flow variation). Possible actions include +, ++, +++, −, −−, and −−−. Here, + (−) means incrementing (decrementing) the base value by one unit. The other possible actions are multiples (twice or three times) of the basic actions + and −.

  • • 

    Number of particles per iteration, #p: This parameter determines the exact number of positions in the binary vector that must be set to one. Base setting: 5; actions: {+, ++, +++, −, −−, −−−}.

  • • 

    Mean and standard deviation: The normal distribution that is sampled in order to set #p positions in the binary vector to one requires a mean—that is, a specific position in the vector—and a standard deviation. Base settings: 12 (mean) and 3 (standard deviation); actions: {+, ++, +++} for the mean and {+, ++, +++, −, −−} for the standard deviation.

  • • 

    (Particle) range: This parameter establishes at which position (to the left and to the right of the mean) the tails of the normal distribution are truncated, resulting in a truncated normal distribution. This implies that all positions outside this range will keep their value of zero. Base setting: 25; actions: {+, ++, +++, −, −−, −−−}.

  • • 

    Variation: This parameter causes a displacement of the particle range towards the left or towards the right, resulting in an asymmetrically truncated normal distribution. Negative values of this parameter indicate a displacement to the left. Base setting: 1, actions: {+, ++, +++}.

The concepts described above are graphically presented in Figure 3. Also note that, for the generation of the bits, the vector is regarded in cyclical form. Performing any of the actions outlined in the list above results in a variation of the type of output vector that is generated, that is, a flow variation. As there are 23 different actions, the flow generation function is able to produce 23 flow variations with respect to the output that is produced by using the base settings.

Figure 3. 

Graphical representation of the parameters used in the context of the flow generation function.

Figure 3. 

Graphical representation of the parameters used in the context of the flow generation function.

In the following, we introduce two different measures related to the output of the flow generation function. At the basis of both measures is the Hamming distance between binary vectors. The following measures are considered:

  • • 
    Rate of change (rc(t, t + 1)): Considering the binary output vector generated at iteration t, and the binary output vector generated at iteration t + 1, the rate of change between iterations t and t + 1 is simply defined as the Hamming distance between and :
    formula
  • • 
    Stability (st(Δ(t, t + i)): While the Hamming distance between two consecutive output vectors provides information about the number of differences in the particles that are generated, it does not provide information about the general stability of the generated particle flow. Therefore, the stability of the flow in a time window Δ(t, t + i)—that is, concerning the iterations between t and t + i—is defined as follows:
    formula

In order to demonstrate the abilities of the flow generation function, we sampled it with the base settings (apart from the standard deviation, which was initially set to zero) for 2000 iterations. Moreover, the standard deviation was increased every 200 iterations until it reached a value of two. The results for the rate of change (between consecutive iterations) and the stability (with a time window size of 100) are shown in Figure 4. In Figure 4a we can observe that the rate of change is zero for the first 200 iterations, which is due to the standard deviation having a value of zero. The same figure shows that, as the standard deviation is increased, the Hamming distance between the generated binary vectors progressively increases as well. Furthermore, Figure 4b shows that in the initial stages of the experiment the stability is at its peak, that is, all generated vectors are very similar. When the standard deviation is increased, the observed stability between consecutive time windows (of size 100) generally decreases.

Figure 4. 

Results demonstrating the characteristics of the flow generation function. (a) Evolution of the rate of change (dH()) between consecutive output vectors generated by the flow generation function over 2000 iterations when the standard deviation is progressively increased from zero to two. (b) Evolution of the stability (st(t, t + 100)) between consecutive time windows of size 100, corresponding to the evolution of the rate of change as shown in (a).

Figure 4. 

Results demonstrating the characteristics of the flow generation function. (a) Evolution of the rate of change (dH()) between consecutive output vectors generated by the flow generation function over 2000 iterations when the standard deviation is progressively increased from zero to two. (b) Evolution of the stability (st(t, t + 100)) between consecutive time windows of size 100, corresponding to the evolution of the rate of change as shown in (a).

For the calibration of the parameter values of the system we used a sequence of 23 flow variations produced by invoking (step by step) the 23 actions outlined in the list above. The initial setting for the flow generation function was always the base setting (also described above). In order to demonstrate the difference between the binary output vectors generated by the 23 flow variations, we generated for each of the 23 settings 100 binary output vectors. The results for the average Hamming distance between these 100 output vectors for each of the 23 settings are shown in Figure 5.

Figure 5. 

The average Hamming distance of 100 binary output vectors as generated by the flow generation function for the 23 different settings used for the system calibration.

Figure 5. 

The average Hamming distance of 100 binary output vectors as generated by the flow generation function for the 23 different settings used for the system calibration.

3.3 Effect of Changing Important System Parameters

In order to demonstrate the importance of calibrating the system parameters, we briefly study two cases. The first one deals with the parameter ϵ (see Section 3.1), which determines the energy loss of agents in the event of a lack of genetic suitability with respect to the particle with which they interact. When the value of this parameter increases, the selection pressure on the population becomes greater, causing an increased sensitivity to flow variations, as shown in Figure 6. This figure also shows that the increase of this parameter value results in a decrease of the population diversity, which is defined as the diversity of the genetic material of the agents of an iteration (see Equation 3). The experiments for Figure 6 were produced with the baseline parameter setting of the system that is outlined in the next section.

Figure 6. 

Evolution of the sensitivity and the population diversity for a changing value of parameter ϵ (see Section 3.1 for a description of the parameter).

Figure 6. 

Evolution of the sensitivity and the population diversity for a changing value of parameter ϵ (see Section 3.1 for a description of the parameter).

The second example concerns the parameter r, whose value determines the minimum energy needed by an agent to reproduce (once again, see Section 3.1 for a description). An increase of the value of this parameter introduces a certain elitism [36] in the population, causing a decrease of population diversity. At the same time the sensitivity increases. Results are shown in Figure 7 and were produced with the baseline parameter settings of the system (see next section).

Figure 7. 

Evolution of the sensitivity and the population diversity for a changing value of parameter r (see Section 3.1 for a description of the parameter).

Figure 7. 

Evolution of the sensitivity and the population diversity for a changing value of parameter r (see Section 3.1 for a description of the parameter).

3.4 Choice of System Parameter Values

First, after initial experiments by hand, we selected a range of presumably suitable values for each of the eight system parameters described in Section 3.1 These parameter value ranges are shown in Table 1. Hereby, the parameter value in the middle of each range (indicated in bold font) corresponds to the baseline setting.

Table 1. 

Parameter value ranges for the eight system parameters. The central value (indicated in bold) within each range corresponds to the baseline value of each parameter.

ParameterConsidered values
ebirth {30, 40, 50, 60, 70} 
einit {50, 75, 100, 125, 150} 
eproc {50, 75, 100, 125, 150} 
r {150, 175, 200, 225, 250} 
vinit {6, 8, 10, 12, 14} 
m {0.4, 0.5, 0.6, 0.7, 0.8} 
ϕ {70, 80, 90, 100, 110} 
ϵ {3.5, 4, 4.5, 5, 5.5} 
ParameterConsidered values
ebirth {30, 40, 50, 60, 70} 
einit {50, 75, 100, 125, 150} 
eproc {50, 75, 100, 125, 150} 
r {150, 175, 200, 225, 250} 
vinit {6, 8, 10, 12, 14} 
m {0.4, 0.5, 0.6, 0.7, 0.8} 
ϕ {70, 80, 90, 100, 110} 
ϵ {3.5, 4, 4.5, 5, 5.5} 

The system was tested with the following parameter settings. First, the baseline setting was applied. Then, one parameter at a time was changed, while the other parameters assumed their baseline values. This results in five different parameter value configurations per parameter. As there are eight parameters in total, this results in 32 parameter value configurations, in addition to the baseline setting. The system was run for each of these 33 parameter value configurations exactly 10 times (for 2000 iterations). Moreover, the flow generation function was used in each run to successively introduce the 23 flow variations as outlined before. For each run we measured the sensitivity of the system (as defined in Equation 2). The average sensitivity values (over 10 runs) per parameter value configuration are presented in the form of barplots in Figure 8. Moreover, this graphic presents the results of the parameter value configuration that is obtained when the best value (according to the tuning experiments) for each parameter is used (label BP). It turns out that this configuration is actually the best one, obtaining an average sensitivity value of around 0.62. Finally, the values of the best parameter configuration are provided in Table 2.

Figure 8. 

Average sensitivity values for the 33 considered parameter value configurations. The last bar (labeled BP) corresponds to the parameter value setting in which the best value is used for each parameter.

Figure 8. 

Average sensitivity values for the 33 considered parameter value configurations. The last bar (labeled BP) corresponds to the parameter value setting in which the best value is used for each parameter.

Table 2. 

Best parameter values found for the FADS in the tuning process.

ParameterBest value
ebirth 50 
einit 150 
eproc 100 
r 250 
vinit 10 
m 0.6 
ϕ 90 
ϵ 5.5 
ParameterBest value
ebirth 50 
einit 150 
eproc 100 
r 250 
vinit 10 
m 0.6 
ϕ 90 
ϵ 5.5 

3.4.1 A Preliminary Application Example concerning Network Security

In the following, we describe the results of a preliminary experiment that we performed with the aim of testing the suitability of the best parameters found for the FADS in the context of network intrusion detection. In particular, we ran the system for 6000 iterations with an initial population of 300 agents equipped with randomly generated genetic code. The input flow of particles was generated on the basis of normal (non-pathological) network traffic generated by Snort in sniffer mode.3 This normal traffic flow was maintained for the first 5570 iterations. In Figure 9 it can be observed that most of the (unfit) agents die during the first 100 iterations, leading to a reduction of the agent population to fewer than 50 individuals. Subsequently, the fit agents proliferate and generate descendants by mutation. Until iteration 5570, we can observe an increase and a stabilization of the agent population (apart from the fractal noise that is intrinsic to ecological systems [24]). In iteration 5570, the normal traffic flow is replaced by particles generated on the basis of a SATAN attack.4 The packets corresponding to the SATAN attack were obtained from the DARPA'98 set of the Lincoln Laboratory at MIT. The evolution of the population size after the attack, as shown in Figure 9, shows a significant fall in the size of the population, that is, a quantifiable deviation in one of the model variables of interest, which will be the basis for using the model as a classifier. Note that a similar decrease in the population size occurs at around iteration 1000. However, at this iteration the system is still in the adaptation phase, eliminating agents that are not adapted and generating offspring of the fittest agents. Moreover, when tstart > 1000 no alert is generated.

Figure 9. 

Evolution of the population size over time.

Figure 9. 

Evolution of the population size over time.

3.5 Population Diversity and Sensitivity

In a second preliminary experiment we aim to show that, with the parameter setting as determined in Section 3.4, after stabilization the agent population is characterized by a low diversity (in terms of the genetic code), which makes the system sensitive to changes in the incoming particle flow. In nature, species lose genetic diversity during evolution due to bottlenecks—such as, for example, the shortage of resources—caused by their environment [39]. This has been shown both by experimental studies [60] and by field studies [18]. These studies show a close relationship between the genetic diversity of a species and its vulnerability to changes in the environment. That is, as genetic diversity diminishes, the species becomes more sensitive to fluctuations in the environment, and its capacity to evolve and adapt to such changes is affected [26]. Although the genetic poverty of a species is a problem that may lead to its extinction, it is actually—as previously mentioned—a desirable characteristic in the context of the FADS developed in this work. This is because an increase in the vulnerability of a species triggers population responses that are more pronounced and easier to detect, and thus more easily allowing the detection of anomalies in the environment.

The data shown in Figure 10 was obtained by running the system starting with 300 agents with random genetic codes (this time using the flow generation function for producing the particle flow) for a total of 33,000 iterations. The standard setting of the flow generation function is used during the first 5000 iterations. It can be observed that, during these 5000 initial iterations, the population diversity (see Section 3 for the exact definition) decreases until a value of approximately 0.15. Then, between iterations 5000 and 23,000, the flow generation function is used for producing variations in the particle flow. It can be seen in Figure 10 that the system reacts with an oscillation of the population diversity. Then, starting from iteration 23,000, the standard settings of the flow generation function are used again for the production of the particle flow, which causes the system to stabilize again, this time at a value of approximately 0.1.

Figure 10. 

Evolution of the population diversity over time.

Figure 10. 

Evolution of the population diversity over time.

4 Adjustment of Classifier Parameters

As previously mentioned in Section 2 and as shown in Algorithm 2, the evolution of the population size (|A(t)|) and the average energy level of the population (e(t)) may be used for the definition of a classifier in order to detect a malicious flow of incoming particles. This classifier has the following five parameters:

  • • 

    ψ: This parameter determines the number of past iterations that are used for the calculation of the average population size SA. See also line 2 of Algorithm 2.

  • • 

    ρ: This parameter determines the number of past iterations that are used for the calculation of the average energy level Se. See also line 3 of Algorithm 2.

  • • 

    dA: A parameter that is used to calculate the bandwidth size (2dA) for the detection of an anomaly. More specifically, an anomaly is detected when |A(t)| surpasses SA + dA or falls below SAdA. See also line 4 of Algorithm 2.

  • • 

    de: A parameter that is used to calculate the bandwidth size (2de) for the detection of an anomaly. More specifically, an anomaly is detected when e(t) surpasses Se + de or falls below Sede. See also line 4 of Algorithm 2.

  • • 

    u: This parameter determines the maximum number of successive iterations over which |A(t)| and/or e(t) may leave their respective bandwidths without triggering an alarm. This concerns the delay mechanism described in Section 2.

First of all, in preliminary experiments we realized that it was only important to detect the cases in which |A(t)| and e(t) fall below their respective lower bounds SAdA and Sede. Therefore, the upper bounds were not considered for the tuning experiments. Moreover, the condition in line 4 of Algorithm 2 can be reduced to the following: an anomaly is detected when |A(t)| < SAdAore(t) < Sede.

During the tuning process, the following values for the above-mentioned parameters were considered:

  • • 

    ψ ∈ {25, 50, 75,…, 900} ∪ {1000, 1200, 1300, 1400, 1500, 4000}

  • • 

    ρ ∈ {25, 50, 75,…, 325} ∪ {400, 600, 800,…, 3200}

  • • 

    dA, de ∈ {0, 1, 2, 4, 8, 16, 32, 64, 128, 256}

  • • 

    u ∈ {0, 1,…, 5} ∪ {10, 15,…, 75} ∪ {100, 200, 300, 600, 1200}

For the purpose of tuning, the system was fed with data packets of 40 different characteristics based on information from the transmission and network layers of the TCP/IP stack. Samples of normal traffic were obtained from the network of the Department of Information Technology of the University of Santiago de Chile, and the security tools Nmap and Nessus were used to generate hostile traffic.5 Specifically, these tools generated probe and denial-of-service (DoS) attacks, respectively. Details concerning the data packets used during the calibration process are provided in Table 3. The best parameter values identified by means of the tuning process are presented in Table 4. These values are used for the final experimental evaluation.

Table 3. 

Data packets used for the tuning of the classifier.

Network traffic source# packets
Normal traffic (University LAN) 37026 
Intrusive traffic (Nessus DoS Scan) 18793 
Intrusive traffic (Nessus Miscellaneous Scan) 85716 
Intrusive traffic (Nmap OS Detection Scan) 14954 
Intrusive traffic (Nmap Syn Scan) 2439 
Intrusive traffic (Nmap Connect Scan) 2073 
Network traffic source# packets
Normal traffic (University LAN) 37026 
Intrusive traffic (Nessus DoS Scan) 18793 
Intrusive traffic (Nessus Miscellaneous Scan) 85716 
Intrusive traffic (Nmap OS Detection Scan) 14954 
Intrusive traffic (Nmap Syn Scan) 2439 
Intrusive traffic (Nmap Connect Scan) 2073 
Table 4. 

Best parameter values found for the classifier in the tuning process.

ParameterBest value
ψ 900 
dA 
ρ 2800 
de 
u 200 
ParameterBest value
ψ 900 
dA 
ρ 2800 
de 
u 200 

5 Experimental Evaluation

The tuning process concerning the proposed FADS, which was described in detail in the previous section, verified our initial research hypothesis, that is, a population of evolutionary agents may become sensitive to environmental changes in a quantifiable manner. Moreover, it even verified that the parameters of the system may be adjusted in order to increase the population sensitivity. In this section, we finally compare the classification results of the FADS with those of other bio-inspired classification techniques that were introduced in the relevant literature for network intrusion detection.

5.1 Experimental Data

The experimental evaluation of the FADS and its competitors was performed on the basis of the widely used DARPA'98 data set. We deliberately did not make use of the KDD'99 data set, which is considered to be of questionable validity and usefulness [49], and which has been shown to be of little use for anomaly detection systems [55]. The DARPA'98 test data set was generated by the Cyber Systems and Technology (CST) group of the Lincoln Laboratory at MIT in the late 1990s under the Defense Advanced Research Projects Agency (DARPA) and the US Air Force Research Laboratory. This set was the first standard data set for the evaluation of intrusion detection systems in computer networks. It was obtained on the basis of a simulation of a local area network (LAN) of a US airbase in which attacks are executed at different workstations. The gross data is organized in terms of nine weeks, separated into seven days per week. We used the entire weeks 1–3; days Mon., Tue., Fri., Sat., Sun. of week 4; days Wed., Thu., Sat., Sun. of week 5; days Thu., Sat., Sun. of week 6; and days Sat., Sun. of week 7 for the training of the competitor algorithms. The rest of the days were used for testing. In this context, remember that the proposed FADS does not require any training period. Finally, note that the DARPA'98 set is freely available at [1].

5.2 Competing Algorithms

Concerning the current state of the art for the DARPA'98 data set, it is almost impossible to single out an individual work. This is because many techniques have been proposed with the aim of improving specific algorithms from the literature. A comprehensive comparison among the latest proposals does not exist. In any case, the following works can be regarded as being part of the state of the art. In [8], the authors presented an approach based on support vector machines that has advantages over existing proposals dealing with such machines. Another highly successful approach is the one proposed in [42], based on a K-nearest-neighbor classifier. Furthermore, neural network approaches can be regarded as being standard approaches for network intrusion detection in the context of the DARPA'98 data set. Well-working examples can be found in [15, 37, 43]. The most successful approach from the field of artificial immune systems is the formal immune network (FIN), proposed in [54]. Finally, the application of dynamical system theory for intrusion detection in the context of the DARPA'98 data set was proposed in [32].

In the following, we describe in more detail the algorithms that were selected for the comparison with the FADS proposed in this work, in addition to the reasons for choosing them:

  • • 

    Artificial neural network (ANN): Because—as mentioned above—neural networks can be regarded as being standard approaches for network intrusion detection, we implemented a classifier based on a multilayer perceptron inspired by [37]. The ANN was trained with backpropagation and the Levenberg-Marquardt algorithm in order to minimize the number of false positives [43]. An initial number of neurons in the hidden layer was calculated according to [27], and then adjusted in a calibration process, leading to the assignment of 16 neurons in the hidden layer.

  • • 

    Formal immune network (FIN): Among the state-of-the-art algorithms for the DARPA'98 data set, FIN [54] is the one that is most closely related to the FADS proposed in this work. The key ingredients of FIN are apoptosis (programmed cell death) and immunization, both of which are controlled by cytokines (messenger proteins). The network is formed by a training process based on the application of a discrete tree transform (DTT) [4] and a singular value decomposition (SVD) [28]. The cytokines operate according to a proximity principle in the FIN space to determine the class of the observed data when mapped into the space by means of the DTT algorithm.

  • • 

    Modified formal immune network (FIN+): Given that both the ANN and FIN are algorithms whose training process works offline, we also wanted to include an algorithm with online learning capabilities into the comparison. For this purpose, we selected a modified (improved) version of FIN named FIN+, which was proposed in [29]. By incorporating the concepts of time and uncertainty, it generates temporary groupings, detecting regions of space with fast growth for the classification as attacks, and thus improving the performance of the original FIN algorithm, particularly in the context of unknown attacks.

5.3 Experimental Design

In order to evaluate the behavior of the four different classifiers in changing scenarios in terms of the detection of known as well as unknown attacks, the following measures were used: (1) the fraction of false positives, henceforth denoted by FP, and (2) a so-called correctness measure as presented in the following equation:
formula
where TP is the number of true positives, TN is the number of true negatives, and n is the total number of analyzed packets. Apart from these two measures, we used the Wilcoxon signed-rank test [61] for the statistical validation of the results.

The tests were divided into two groups: known attacks and unknown attacks. Tests concerning known attacks involved the training of the classifiers ANN, FIN, and FIN+, including attacks to which they were later exposed during testing. In this context, it is important to remember that a FADS does not require a prior training process. Instead, before starting with the attacks, the FADS was exposed to 3000 packets of normal traffic, which is a number empirically shown to be sufficient for the stabilization of the agent population. Note that the system is always under adaptation, not only at the start. It other words, the system is continuously subject to unsupervised adaptation. Obviously, at the beginning the system is vulnerable, but the time required for about 3000 packets is negligible in a data network where packet counts are measured in millions of packets per second (Mpps). Finally, note that tests concerning unknown attacks are based on attacks not included in the training sets of ANN, FIN, and FIN+. Table 5 shows the total numbers of normal packets and attack packets used in the tests, and the numbers of attacks per week.6

Table 5. 

Total number of packets used for testing.

Week# of normal packets# of attack packets# of attacks
1,377,730 632 
1,548,789 62,407 
876,942 315,294 13 
1,114,619 5,220,620 
1,302,544 541,326 
1,627,996 879,216 12 
1,660,560 264,181 11 
Test 1 2,032,815 674,120 17 
Test 2 2,445,865 1,229,280 44 
Week# of normal packets# of attack packets# of attacks
1,377,730 632 
1,548,789 62,407 
876,942 315,294 13 
1,114,619 5,220,620 
1,302,544 541,326 
1,627,996 879,216 12 
1,660,560 264,181 11 
Test 1 2,032,815 674,120 17 
Test 2 2,445,865 1,229,280 44 

6 Results

The numerical results are presented in a summarized form in Table 6. In particular, the table provides the average fraction of false positives (FPs) and the value of the correctness measure (as introduced in Equation 6) obtained by the four different classifiers for tests with known and unknown attacks. It may be observed that, in the context of known attacks, the results of the FADS (15.56%) are characterized by fewer false positives than those of the other classifiers. The second-best algorithm is the modified FIN version by Iturbe (FIN+) [29], which obtains 25.38% false positives. The results of the Wilcoxon signed-rank test indicate that the two algorithms are different, with a p-value smaller than 0.05. Moreover, FADS also outperforms its competitors in terms of the correctness measure. However, in this case the differences between the FADS and the second-best approach (FIN+) are not statistically significant, that is, we obtained a p-value greater than 0.05.

Table 6. 

Summarized results. Detailed numerical results are shown in Table 7.

ClassifierKnown attacksUnknown attacks
Fraction of FPsCorrectnessFraction of FPsCorrectness
FADS 0.1556 0.7848 0.1556 0.7848 
FIN 0.2600 0.7017 0.2923 0.6756 
FIN+ 0.2538 0.6728 0.2762 0.6798 
ANN 0.3097 0.6670 0.3477 0.6284 
ClassifierKnown attacksUnknown attacks
Fraction of FPsCorrectnessFraction of FPsCorrectness
FADS 0.1556 0.7848 0.1556 0.7848 
FIN 0.2600 0.7017 0.2923 0.6756 
FIN+ 0.2538 0.6728 0.2762 0.6798 
ANN 0.3097 0.6670 0.3477 0.6284 

With regard to unknown attacks, we must mention, first of all, that the results of the FADS are exactly the same as in the case of known attacks. This is because the FADS does not require an initial training phase and therefore all attacks are unknown. In contrast to the FADS, the results of the remaining classifiers naturally become worse when moving to unknown attacks. In fact, the FADS outperforms the competitors now in a statistically significant way, both in terms of the fraction of false positives and in terms of the correctness measure. Detailed numerical results (for each day that was considered for testing) are provided in Table 7. The best results for each day are entered in bold. The last row in the table provides once more the summarized results that were previously presented in Table 6.

Table 7. 

Detailed numerical results for ANN, FIN, FIN+, and FADS concerning both the fraction of false positives (FP) and correctness in the presence of known and unknown attacks.

WeekDayKnown attacksUnknown attacks
Fraction of FPsAccuracyFraction of FPsAccuracy
ANNFINFIN+FADSANNFINFIN+FADSANNFINFIN+FADSANNFINFIN+FADS
Wednesday 0.6165 0.3996 0.3248 0.1456 0.3835 0.6000 0.6748 0.8528 0.6185 0.4008 0.3258 0.1456 0.3815 0.5992 0.6742 0.8528 
Thursday 0.5851 0.5590 0.5064 0.2071 0.3908 0.4173 0.4712 0.7401 0.6230 0.5951 0.5339 0.2071 0.3770 0.4049 0.4661 0.7401 
Monday 0.5286 0.4508 0.3839 0.2121 0.4612 0.5426 0.6096 0.7825 0.5538 0.4722 0.3897 0.2121 0.4462 0.5278 0.6103 0.7825 
Tuesday 0.4332 0.2619 0.7682 0.1430 0.5668 0.7381 0.2318 0.8550 0.4342 0.2625 0.7699 0.1430 0.5658 0.7375 0.2301 0.8550 
Friday 0.3053 0.2277 0.2032 0.2179 0.6938 0.7521 0.7749 0.7273 0.5906 0.4403 0.3931 0.2179 0.4094 0.5597 0.6069 0.7273 
Monday 0.5220 0.4343 0.4131 0.3665 0.4773 0.5101 0.3182 0.5785 0.5682 0.4727 0.3893 0.3665 0.4318 0.5273 0.6107 0.5785 
Tuesday 0.5019 0.4098 0.3514 0.2571 0.4981 0.5876 0.9576 0.7016 0.5351 0.4368 0.3718 0.2571 0.4649 0.5632 0.6282 0.7016 
Wednesday 0.2332 0.2009 0.2054 0.2359 0.7581 0.7758 0.9406 0.7251 0.3794 0.3268 0.2962 0.2359 0.6077 0.6604 0.6910 0.7251 
Friday 0.3049 0.2504 0.2213 0.1672 0.6951 0.7443 0.0965 0.8246 0.3067 0.2518 0.2174 0.1672 0.6933 0.7482 0.7826 0.8246 
Monday 0.3224 0.2562 0.2308 0.1157 0.6776 0.7438 0.1111 0.8842 0.3224 0.2562 0.2307 0.1157 0.6776 0.7438 0.7693 0.8842 
Tuesday 0.2893 0.2830 0.2461 0.1389 0.7107 0.7167 0.9303 0.8596 0.2906 0.2843 0.2469 0.1389 0.7094 0.7157 0.7531 0.8596 
Wednesday 0.2927 0.2481 0.2114 0.1456 0.7073 0.7517 0.9296 0.8524 0.2934 0.2486 0.2117 0.1456 0.7066 0.7514 0.7883 0.8524 
Thursday 0.3063 0.2553 0.2172 0.1288 0.6899 0.7425 0.8317 0.5908 0.3122 0.2601 0.2181 0.1288 0.6878 0.7399 0.7819 0.5908 
Friday 0.2269 0.1638 0.1540 0.1151 0.7655 0.8166 0.9412 0.7998 0.3367 0.2431 0.1981 0.1151 0.6532 0.7469 0.7919 0.7998 
Test 1 Monday 0.1769 0.1992 0.1650 0.0706 0.8231 0.8006 0.8348 0.9282 0.1769 0.1993 0.1651 0.0706 0.8231 0.8007 0.8349 0.9282 
Tuesday 0.1565 0.1662 0.1401 0.1208 0.8269 0.8030 0.8280 0.7980 0.2294 0.2434 0.2007 0.1208 0.7504 0.7323 0.7750 0.7980 
Wednesday 0.2740 0.2322 0.2179 0.1024 0.7244 0.7675 0.7817 0.7344 0.2744 0.2324 0.2182 0.1024 0.7243 0.7674 0.7816 0.7344 
Thursday 0.2709 0.2359 0.2213 0.0550 0.6628 0.7040 0.7187 0.9450 0.2710 0.2361 0.2214 0.0550 0.6627 0.7040 0.7187 0.9450 
Friday 0.1220 0.1284 0.1088 0.1631 0.7916 0.7683 0.7879 0.7399 0.1844 0.1940 0.1645 0.1631 0.6862 0.6723 0.7019 0.7399 
Test 2 Monday 0.1691 0.1626 0.1666 0.1173 0.7287 0.6966 0.7707 0.8234 0.1713 0.1647 0.1688 0.1173 0.7252 0.7023 0.7773 0.8234 
Tuesday 0.1995 0.1549 0.1487 0.1092 0.7758 0.7609 0.6749 0.8044 0.2081 0.1616 0.1548 0.1092 0.7720 0.7667 0.6762 0.8044 
Wednesday 0.1418 0.1312 0.1143 0.1142 0.7840 0.7203 0.6487 0.7671 0.1475 0.1364 0.1189 0.1142 0.7755 0.7138 0.6392 0.7671 
Thursday 0.1675 0.1958 0.1426 0.1256 0.7614 0.7154 0.7086 0.7520 0.2165 0.2532 0.1844 0.1256 0.6924 0.6489 0.6410 0.7520 
Friday 0.2874 0.2324 0.2289 0.1604 0.6532 0.6647 0.5737 0.7681 0.3003 0.2428 0.2392 0.1604 0.6581 0.6798 0.5843 0.7681 
Average 0.3097 0.2600 0.2538 0.1556 0.6670 0.7017 0.6728 0.7848 0.3477 0.2923 0.2762 0.1556 0.6284 0.6756 0.6798 0.7848 
WeekDayKnown attacksUnknown attacks
Fraction of FPsAccuracyFraction of FPsAccuracy
ANNFINFIN+FADSANNFINFIN+FADSANNFINFIN+FADSANNFINFIN+FADS
Wednesday 0.6165 0.3996 0.3248 0.1456 0.3835 0.6000 0.6748 0.8528 0.6185 0.4008 0.3258 0.1456 0.3815 0.5992 0.6742 0.8528 
Thursday 0.5851 0.5590 0.5064 0.2071 0.3908 0.4173 0.4712 0.7401 0.6230 0.5951 0.5339 0.2071 0.3770 0.4049 0.4661 0.7401 
Monday 0.5286 0.4508 0.3839 0.2121 0.4612 0.5426 0.6096 0.7825 0.5538 0.4722 0.3897 0.2121 0.4462 0.5278 0.6103 0.7825 
Tuesday 0.4332 0.2619 0.7682 0.1430 0.5668 0.7381 0.2318 0.8550 0.4342 0.2625 0.7699 0.1430 0.5658 0.7375 0.2301 0.8550 
Friday 0.3053 0.2277 0.2032 0.2179 0.6938 0.7521 0.7749 0.7273 0.5906 0.4403 0.3931 0.2179 0.4094 0.5597 0.6069 0.7273 
Monday 0.5220 0.4343 0.4131 0.3665 0.4773 0.5101 0.3182 0.5785 0.5682 0.4727 0.3893 0.3665 0.4318 0.5273 0.6107 0.5785 
Tuesday 0.5019 0.4098 0.3514 0.2571 0.4981 0.5876 0.9576 0.7016 0.5351 0.4368 0.3718 0.2571 0.4649 0.5632 0.6282 0.7016 
Wednesday 0.2332 0.2009 0.2054 0.2359 0.7581 0.7758 0.9406 0.7251 0.3794 0.3268 0.2962 0.2359 0.6077 0.6604 0.6910 0.7251 
Friday 0.3049 0.2504 0.2213 0.1672 0.6951 0.7443 0.0965 0.8246 0.3067 0.2518 0.2174 0.1672 0.6933 0.7482 0.7826 0.8246 
Monday 0.3224 0.2562 0.2308 0.1157 0.6776 0.7438 0.1111 0.8842 0.3224 0.2562 0.2307 0.1157 0.6776 0.7438 0.7693 0.8842 
Tuesday 0.2893 0.2830 0.2461 0.1389 0.7107 0.7167 0.9303 0.8596 0.2906 0.2843 0.2469 0.1389 0.7094 0.7157 0.7531 0.8596 
Wednesday 0.2927 0.2481 0.2114 0.1456 0.7073 0.7517 0.9296 0.8524 0.2934 0.2486 0.2117 0.1456 0.7066 0.7514 0.7883 0.8524 
Thursday 0.3063 0.2553 0.2172 0.1288 0.6899 0.7425 0.8317 0.5908 0.3122 0.2601 0.2181 0.1288 0.6878 0.7399 0.7819 0.5908 
Friday 0.2269 0.1638 0.1540 0.1151 0.7655 0.8166 0.9412 0.7998 0.3367 0.2431 0.1981 0.1151 0.6532 0.7469 0.7919 0.7998 
Test 1 Monday 0.1769 0.1992 0.1650 0.0706 0.8231 0.8006 0.8348 0.9282 0.1769 0.1993 0.1651 0.0706 0.8231 0.8007 0.8349 0.9282 
Tuesday 0.1565 0.1662 0.1401 0.1208 0.8269 0.8030 0.8280 0.7980 0.2294 0.2434 0.2007 0.1208 0.7504 0.7323 0.7750 0.7980 
Wednesday 0.2740 0.2322 0.2179 0.1024 0.7244 0.7675 0.7817 0.7344 0.2744 0.2324 0.2182 0.1024 0.7243 0.7674 0.7816 0.7344 
Thursday 0.2709 0.2359 0.2213 0.0550 0.6628 0.7040 0.7187 0.9450 0.2710 0.2361 0.2214 0.0550 0.6627 0.7040 0.7187 0.9450 
Friday 0.1220 0.1284 0.1088 0.1631 0.7916 0.7683 0.7879 0.7399 0.1844 0.1940 0.1645 0.1631 0.6862 0.6723 0.7019 0.7399 
Test 2 Monday 0.1691 0.1626 0.1666 0.1173 0.7287 0.6966 0.7707 0.8234 0.1713 0.1647 0.1688 0.1173 0.7252 0.7023 0.7773 0.8234 
Tuesday 0.1995 0.1549 0.1487 0.1092 0.7758 0.7609 0.6749 0.8044 0.2081 0.1616 0.1548 0.1092 0.7720 0.7667 0.6762 0.8044 
Wednesday 0.1418 0.1312 0.1143 0.1142 0.7840 0.7203 0.6487 0.7671 0.1475 0.1364 0.1189 0.1142 0.7755 0.7138 0.6392 0.7671 
Thursday 0.1675 0.1958 0.1426 0.1256 0.7614 0.7154 0.7086 0.7520 0.2165 0.2532 0.1844 0.1256 0.6924 0.6489 0.6410 0.7520 
Friday 0.2874 0.2324 0.2289 0.1604 0.6532 0.6647 0.5737 0.7681 0.3003 0.2428 0.2392 0.1604 0.6581 0.6798 0.5843 0.7681 
Average 0.3097 0.2600 0.2538 0.1556 0.6670 0.7017 0.6728 0.7848 0.3477 0.2923 0.2762 0.1556 0.6284 0.6756 0.6798 0.7848 

Finally, we studied the reasons for which the FADS was not able to detect some of the attacks. In the course of these investigations we realized that none of these attacks was long enough (in terms of the number of system iterations) to cause the launching of an alert message, because the value of the variable u—responsible for the delay of alert messages—was set during tuning to 200. Thus, any attack that was not part of a packet trace (with normal background traffic)—that is, longer than 200—failed to trigger an alarm for the classifier. This may be related to the fact that tuning was performed with probe and DoS attacks of moderate length. It should be noted that failure to detect these short attacks in the tests did not significantly affect the results, owing to the low number of short attacks in the test data set of DARPA'98.

7 Conclusions

In this paper, we have presented a novel model for anomaly detection inspired by an ecological view of biological immune systems. The basic idea is that the evolution of a population of agents, while they adapt to changes in the observed data, allows one to discover anomalous behaviors in the system. To demonstrate its validity, the proposed model has been applied to a network intrusion detection problem. In this particular application there are mainly two characteristics that distinguish our model. Firstly, it is able to discover previous unseen attacks, and secondly (in contrast to most of the proposals in the literature), it does not require any previous training. A set of experiments have been carried out comparing our proposal with state-of-the-art algorithms on common benchmarks. The results confirm the validity of our approach. Furthermore, its successes in discovering unseen attacks are outstanding.

Our proposal can be seen as an initial step toward the design of models based on this biological tendency of the ecological view to immune systems. As a pioneer approach, it may certainly be subject to a series of improvements. A first improvement would be the reduction of the number of parameters of the system. However this should be done while keeping as close as possible to the source of inspiration. The calibration process is a second opportunity for improvement. Currently this process has been carried out with a trial-and-error heuristic algorithm. It is necessary to conduct a study to clearly define the design and appropriate principles of these processes.

Acknowledgments

C. Blum was supported by project TIN2012-37930 of the Spanish Government. In addition, support is acknowledged from IKERBASQUE (Basque Foundation for Science). J. A. Lozano was partially supported by the IT609-13 program (Basque Government) and TIN2013-41272-P (Spanish Ministry of Science and Innovation).

Notes

1 

The von Neumann neighborhood of a grid cell consists of the four cells orthogonally surrounding the cell.

2 

Remember that the parameter u is used in the context of the time-delayed alert mechanism implemented by Algorithm 2.

3 

Snort is a free and open-source network intrusion prevention system and network intrusion detection system; see www.snort.org.

4 

SATAN is the abbreviation for Security Administrator Tool for Analyzing Networks.

5 

Nmap (http://nmap.org) is an open source tool for network exploration and security auditing. It was designed to rapidly scan large networks, but it also works fine in the context of single hosts. Nessus (http://www.tenable.com/products/nessus) is a vulnerability and configuration assessment product.

6 

Remember that the DARPA'98 set is divided into days and weeks.

References

1
DARPA'98 Data Set
. .
2
Aickelin
,
U.
,
Greensmith
,
J.
, &
Twycross
,
J.
(
2004
).
Immune system approaches to intrusion detection—a review
. In
G.
Nicosia
,
V.
Cutello
,
P. J.
Bentley
, &
J.
Timmis
(Eds.),
Proceedings of Artificial Immune Systems
(pp.
316
329
).
Cham, Heidelberg, New York, Dordrecht, London
:
Springer
.
3
Amor
,
N. B.
,
Benferhat
,
S.
, &
Elouedi
,
Z.
(
2004
).
Naive Bayes vs decision trees in intrusion detection systems
. In
Proceedings of SAC 2004—ACM Symposium on Applied Computing
(pp.
420
424
).
New York
:
ACM Press
.
4
Atreas
,
N. D.
,
Karanikas
,
C. G.
, &
Tarakanov
,
A. O.
(
2003
).
Signal processing by an immune type tree transform
. In
J.
Timmis
,
P. J.
Bentley
, &
E.
Hart
(Eds.),
Proceedings of Artificial Immune Systems
(pp.
111
119
).
Cham, Heidelberg, New York, Dordrecht, London
:
Springer
.
5
Bersini
,
H.
(
2002
).
Self-assertion versus self-recognition: A tribute to Francisco Varela
. In
J.
Timmis
&
P. J.
Bentley
,
(Eds.),
Proceedings of ICARIS 2002—1st International Conference on Artificial Immune Systems
(pp.
107
112
).
Canterbury
:
University of Kent at Canterbury Printing Unit
.
6
Bersini
,
H.
, &
Varela
,
F. J.
(
1991
).
Hints for adaptive problem solving gleaned from immune network
. In
H.-P.
Schwefel
&
R.
Männer
(Eds.),
Proceedings of PPSN 1991—Parallel Problem Solving from Nature
(pp.
343
354
).
Cham, Heidelberg, New York, Dordrecht, London
:
Springer
.
7
Catania
,
C.
, &
Garcia
,
C.
(
2012
).
Automatic network intrusion detection: Current techniques and open issues
.
Computers & Electrical Engineering
,
38
(
5
),
1062
1072
.
8
Catania
,
C. A.
,
Bromberg
,
F.
, &
García Garino
,
C.
(
2012
).
An autonomous labeling approach to support vector machines algorithms for network traffic anomaly detection
.
Expert Systems with Applications
,
39
(
2
),
1822
1829
.
9
Cho
,
S. B.
(
2002
).
Incorporating soft computing techniques into a probabilistic intrusion detection system
.
IEEE Transactions on Systems, Man, and Cybernetics Part C—Applications and Reviews
,
32
(
2
),
154
160
.
10
Ciccazzo
,
A.
,
Conca
,
P.
,
Nicosia
,
G.
, &
Stracquadanio
,
G.
(
2008
).
An advanced clonal selection algorithm with ad-hoc network-based hypermutation operators for synthesis of topology and sizing of analog electrical circuits
. In
P. J.
Bentley
,
D.
Lee
, &
S.
Jung
(Eds.),
Proceedings of Artificial Immune Systems
(pp.
60
70
).
Cham, Heidelberg, New York, Dordrecht, London
:
Springer
.
11
Dasgupta
,
D.
(
1998
).
Artificial immune systems and their applications
.
Berlin, Heidelberg
:
Springer
.
12
De Castro
,
L.
, &
Von Zuben
,
F.
(
2000
).
The clonal selection algorithm with engineering applications
. In
Proceedings of the GECCO Workshop on Artificial Immune Systems and Their Applications
(pp.
36
37
).
New York
:
ACM Press
.
13
De Castro
,
L. N.
, &
Timmis
,
J.
(
2002
).
Artificial immune systems: A new computational intelligence approach
.
Cham, Heidelberg, New York, Dordrecht, London
:
Springer
.
14
De Castro
,
L. N.
, &
Zuben
,
F. J. V.
(
2001
).
aiNET: An artificial immune network for data analysis
. In
H. A.
Abbass
,
R. A.
Sarker
, &
C. S.
Newton
(Eds.),
Data mining: A heuristic approach
(pp.
231
259
).
Hershey
:
Idea Group Publishing
.
15
Elfeshawy
,
N. A.
, &
Faragallah
,
O. S.
(
2013
).
Divided two-part adaptive intrusion detection system
.
Wireless Networks
,
19
(
3
),
301
321
.
16
Estevez-Tapiador
,
J. M.
,
Garcia-Teodoro
,
P.
, &
Diaz-Verdejo
,
J. E.
(
2004
).
Anomaly detection methods in wired networks: A survey and taxonomy
.
Computer Communications
,
27
(
16
),
1569
1584
.
17
Farmer
,
J. D.
,
Packard
,
N. H.
, &
Perelson
,
A. S.
(
1986
).
The immune system, adaptation, and machine learning
.
Physica D: Nonlinear Phenomena
,
22
(
1
),
187
204
.
18
Fauvelot
,
C.
,
Bernardi
,
G.
, &
Planes
,
S.
(
2003
).
Reductions in the mitochondrial DNA diversity of coral reef fish provide evidence of population bottlenecks resulting from Holocene sea-level change
.
Evolution
,
57
(
7
),
1571
1583
.
19
Forrest
,
S.
,
Perelson
,
A. S.
,
Allen
,
L.
, &
Cherukuri
,
R.
(
1994
).
Self-nonself discrimination in a computer
. In
Proceedings of the 2012 IEEE Symposium on Security and Privacy
(pp.
202
212
).
Washington
:
IEEE Computer Society Press
.
20
Garcia-Teodoro
,
P.
,
Diaz-Verdejo
,
J.
,
Macia-Fernandez
,
G.
, &
Vazquez
,
E.
(
2009
).
Anomaly-based network intrusion detection: Techniques, systems and challenges
.
Computers & Security
,
28
(
1–2
),
18
28
.
21
Greensmith
,
J.
,
Aickelin
,
U.
, &
Twycross
,
J.
(
2010
).
Detecting danger: Applying a novel immunological concept to intrusion detection systems
.
ArXiv preprint arXiv:1002.0696
.
22
Greitzer
,
F. L.
,
Moore
,
A. P.
,
Cappelli
,
D. M.
,
Andrews
,
D. H.
,
Carroll
,
L. A.
, &
Hull
,
T. D.
(
2008
).
Combating the insider cyber threat
.
IEEE Security & Privacy
,
6
(
1
),
61
64
.
23
Halavati
,
R.
,
Shouraki
,
S. B.
,
Heravi
,
M. J.
, &
Jashmi
,
B. J.
(
2007
).
An artificial immune system with partially specified antibodies
. In
Proceedings of the 9th Annual Conference on Genetic and Evolutionary Computation
(pp.
57
62
).
New York
:
ACM Press
.
24
Halley
,
J. M.
(
1996
).
Ecology, evolution and 1/f-noise
.
Trends in Ecology & Evolution
,
11
(
1
),
33
37
.
25
Harmer
,
P.
,
Williams
,
P.
,
Gunsch
,
G.
, &
Lamont
,
G.
(
2002
).
An artificial immune system architecture for computer security applications
.
IEEE Transactions on Evolutionary Computation
,
6
(
3
),
252
280
.
26
Hauser
,
L.
,
Adcock
,
G.
,
Smith
,
P.
,
Ramirez
,
J.
, &
Carvalho
,
G.
(
2002
).
Loss of microsatellite diversity and low effective population size in an overexploited population of New Zealand snapper (Pagrus auratus)
.
Proceedings of the National Academy of Sciences
,
99
(
18
),
11742
11747
.
27
Haykin
,
S.
(
2008
).
Neural networks and learning machines
(3rd ed.).
Upper Saddle River, NJ
:
Prentice Hall
.
28
Horn
,
R.
, &
Johnson
,
C.
(
1986
).
Matrix analysis
.
Cambridge, UK
:
Cambridge University Press
.
29
Iturbe
,
J. A.
(
2010
).
Adaptation of algorithms inspired by artificial immune systems for dynamic intrusion detection in data networks
(in Spanish). Unpublished master's thesis, Universidad de Santiago de Chile, Chile
.
30
Jeffrey
,
D. W.
, &
Madden
,
B.
(
1991
).
Bioindicators and environmental management
.
London
:
Academic Press
.
31
Ji
,
Z.
, &
Dasgupta
,
D.
(
2007
).
Revisiting negative selection algorithms
.
Evolutionary Computation
,
15
(
2
),
223
251
.
32
Kanaskar
,
N.
,
Seker
,
R.
,
Jiang
,
B.
, &
Phoha
,
V. V.
(
2012
).
Dynamical system theory for the detection of anomalous behavior in computer programs
.
IEEE Transactions on Systems, Man, and Cybernetics, Part C: Applications and Reviews
,
42
(
6
),
1579
1589
.
33
Kephart
,
J. O.
(
1994
).
A biologically inspired immune system for computers
. In
R. A.
Brooks
&
P.
Maes
(Eds.),
Proceedings of Artificial Life IV—Fourth International Workshop on the Synthesis and Simulation of Living Systems
(pp.
130
139
).
Cambridge, MA
:
Massachusetts Institute of Technology
.
34
Kim
,
J.
, &
Bentley
,
P. J.
(
2001
).
Towards an artificial immune system for network intrusion detection: An investigation of clonal selection with a negative selection operator
. In
Proceedings of the 2001 Congress on Evolutionary Computation
,
Vol. 2
(pp.
1244
1252
).
Washington
:
IEEE Press
.
35
Kolias
,
C.
,
Kambourakis
,
G.
, &
Maragoudakis
,
M.
(
2011
).
Swarm intelligence in intrusion detection: A survey
.
Computers & Security
,
30
(
8
),
625
642
.
36
Koza
,
J. R.
(
1992
).
Genetic programming: On the programming of computers by means of natural selection
.
Cambridge, MA
:
MIT Press
.
37
Kukielka
,
P.
, &
Kotulski
,
Z.
(
2008
).
Analysis of different architectures of neural networks for application in intrusion detection systems
. In
Proceedings of IMCSIT 2008—International Multiconference on Computer Science and Information Technology
(pp.
807
811
).
Washington
:
IEEE Press
.
38
Kurlej
,
B.
, &
Wozniak
,
M.
(
2012
).
Active learning approach to concept drift problem
.
Logic Journal of the IGPL
,
20
(
3
),
550
559
.
39
Lambert
,
D.
,
King
,
T.
,
Shepherd
,
L.
,
Livingston
,
A.
,
Anderson
,
S.
, &
Craig
,
J.
(
2005
).
Serial population bottlenecks and genetic variation: Translocated populations of the New Zealand Saddleback (Philesturnus carunculatus rufusater)
.
Conservation Genetics
,
6
(
1
),
1
14
.
40
Li
,
Y.
,
Xia
,
J.
,
Zhang
,
S.
,
Yan
,
J.
,
Ai
,
X.
, &
Dai
,
K.
(
2012
).
An efficient intrusion detection system based on support vector machines and gradually feature removal method
.
Expert Systems with Applications
,
39
(
1
),
424
430
.
41
Liao
,
H. J.
,
Lin
,
C. H.
,
Lin
,
Y. C.
, &
Tung
,
K. Y.
(
2013
).
Intrusion detection system: A comprehensive review
.
Journal of Network and Computer Applications
,
36
(
1
),
16
24
.
42
Liao
,
Y.
, &
Vemuri
,
V. R.
(
2002
).
Use of K-nearest neighbor classifier for intrusion detection
.
Computers & Security
,
21
(
5
),
439
448
.
43
Linda
,
O.
,
Vollmer
,
T.
, &
Manic
,
M.
(
2009
).
Neural network based intrusion detection system for critical infrastructures
. In
Proceedings of IJCNN 2009—International Joint Conference on Neural Networks
(pp.
1827
1834
).
Washington
:
IEEE Press
.
44
Lu
,
W.
, &
Traore
,
I.
(
2004
).
Detecting new forms of network intrusion using genetic programming
.
Computational Intelligence
,
20
(
3
),
475
494
.
45
Matzinger
,
P.
(
2002
).
The danger model: A renewed sense of self
.
Science
,
296
(
5566
),
301
305
.
46
Mostardinha
,
P.
,
Faria
,
B. F.
,
Zúquete
,
A.
, &
Vistulo de Abreu
,
F.
(
2012
).
A negative selection approach to intrusion detection
. In
C. A.
Coello-Coello
,
J.
Greensmith
,
N.
Krasnogor
,
P.
Liò
,
G.
Nicosia
, &
M.
Pavone
(Eds.),
Proceedings of Artificial Immune Systems
(pp.
178
190
).
Cham, Heidelberg, New York, Dordrecht, London
:
Springer
.
47
Mukherjee
,
B.
,
Heberlein
,
L.
, &
Levitt
,
K.
(
1994
).
Network intrusion detection
.
IEEE Network
,
8
(
3
),
26
41
.
48
Newman
,
D.
,
Snyder
,
J.
, &
Thayer
,
R.
(
2002
).
Crying wolf: False alarms hide attacks
Network World
, .
49
Olusola
,
A.
,
Oladele
,
A.
, &
Abosede
,
D.
(
2010
).
Analysis of KDD'99 intrusion detection dataset for selection of relevance features
. In
S. I.
Ao
,
C.
Douglas
,
W. S.
Grundfest
, &
J.
Burgstone
(Eds.),
Proceedings of the 2010 World Congress on Engineering and Computer Science
(pp.
162
168
).
Hong Kong
:
Newswood Limited
.
50
Ramadas
,
M.
,
Ostermann
,
S.
, &
Tjaden
,
B.
(
2003
).
Detecting anomalous network traffic with self-organizing maps
. In
G.
Vigna
,
E.
Jonsson
, &
C.
Kruegel
(Eds.),
Proceedings of RAID 2003—6th International Symposium on Recent Advances in Intrusion Detection
(pp.
36
54
).
Cham, Heidelberg, New York, Dordrecht, London
:
Springer
.
51
Sklar
,
E.
(
2007
).
NetLogo, a multi-agent simulation environment
.
Artificial Life
,
13
(
3
),
303
311
.
52
Stopel
,
D.
,
Moskovitch
,
R.
,
Boger
,
Z.
,
Shahar
,
Y.
, &
Elovici
,
Y.
(
2009
).
Using artificial neural networks to detect unknown computer worms
.
Neural Computing and Applications
,
18
(
7
),
663
674
.
53
Suzuki
,
H.
,
Ono
,
N.
, &
Yuta
,
K.
(
2003
).
Several necessary conditions for the evolution of complex forms of life in an artificial environment
.
Artificial Life
,
9
(
2
),
153
174
.
54
Tarakanov
,
A.
(
2008
).
Immunocomputing for intelligent intrusion detection
.
IEEE Computational Intelligence Magazine
,
3
(
2
),
22
30
.
55
Tavallaee
,
M.
,
Bagheri
,
E.
,
Wei
,
L.
, &
Ghorbani
,
A. A.
(
2009
).
A detailed analysis of the KDD CUP 99 data set
. In
Proceedings of CISDA 2009—IEEE Symposium on Computational Intelligence for Security and Defense Applications
(pp.
1
6
).
Washington
:
IEEE Press
.
56
Tsai
,
C. F.
,
Hsu
,
Y. F.
,
Lin
,
C. Y.
, &
Lin
,
W. Y.
(
2009
).
Intrusion detection by machine learning: A review
.
Expert Systems with Applications
,
36
(
10
),
11994
12000
.
57
Tsang
,
C. H.
,
Kwong
,
S.
, &
Wang
,
H.
(
2007
).
Genetic-fuzzy rule mining approach and evaluation of feature selection techniques for anomaly intrusion detection
.
Pattern Recognition
,
40
(
9
),
2373
2391
.
58
Varela
,
F.
(
2000
).
El fenómeno de la vida
.
Santiago, Chile
:
Dolmen Editorial
.
59
Vaz
,
N. M.
, &
Varela
,
F. J.
(
1978
).
Self and non-sense: An organism-centered approach to immunology
.
Medical Hypotheses
,
4
(
3
),
231
267
.
60
Wahl
,
L. M.
,
Gerrish
,
P. J.
, &
Saika-Voivod
,
I.
(
2002
).
Evaluating the impact of population bottlenecks in experimental evolution
.
Genetics
,
162
(
2
),
961
971
.
61
Wilcoxon
,
F.
(
1945
).
Individual comparison by ranking methods
.
Biometrics Bulletin
,
1
(
6
),
80
83
.
62
Wu
,
S.
, &
Banzhaf
,
W.
(
2010
).
The use of computational intelligence in intrusion detection systems: A review
.
Applied Soft Computing
,
10
(
1
),
1
35
.
63
Ye
,
N.
,
Emran
,
S. M.
,
Chen
,
Q.
, &
Vilbert
,
S.
(
2002
).
Multivariate statistical analysis of audit trails for host-based intrusion detection
.
IEEE Transactions on Computers
,
51
(
7
),
810
820
.

Author notes

Contact author.

∗∗

Department of Computer Science and Artificial Intelligence, University of the Basque Country UPV/EHU, Paseo Manuel Lardizabal 1, 20018 San Sebastian, Spain. E-mail: christian.blum@ehu.es (C.B.); ja.lozano@ehu.es (J.A.L.)

IKERBASQUE, Basque Foundation for Science, Bilbao, Spain.

Escuela de Informática, Universidad Santo Tomás, Avenida Prat 855, Concepcián, Chile. E-mail: ppinacho@santotomas.cl