Bogusław Cyganek
email: cyganek@agh.edu.pl

Embedding of the Extended Euclidean Distance into Pattern Recognition with Higher-Order Singular Value Decomposition of Prototype Tensors

Keywords: Pattern classification, ensemble of classifiers, Euclidean Distance, IMED, HOSVD

The paper presents architecture and properties of the ensemble of the classifiers operating in the tensor orthogonal spaces obtained with the Higher-Order Singular Value Decomposition of prototype tensors. In this paper two modifications to this architecture are proposed. The first one consists in embedding of the Extended Euclidean Distance metric which accounts for the spatial relationship of pixels in the input images and allows robustness to small geometrical perturbations of the patterns. The second improvement consists in application of the weighted majority voting for combination of the responses of the classifiers in the ensemble. The experimental results show that the proposed improvements increase overall accuracy of the ensemble.

Introduction

This paper is an extension of our previous work on development of the image classification with the ensemble of tensor based classifiers [START_REF] Cyganek | Ensemble of Tensor Classifiers Based on the Higher-Order Singular Value Decomposition[END_REF]. The method showed to be very robust in terms of accuracy and execution time, since many existing methods do not account for the multi-dimensionality of the classified data [START_REF] Szeliski | Computer Vision. Algorithms and Applications[END_REF] [START_REF] Turk | Eigenfaces for recognition[END_REF] [START_REF] Vasilescu | Multilinear analysis of image ensembles: TensorFaces[END_REF].

Processing and classification of the multi-factor dependent data can be addressed with help of methods operating with tensors and their decompositions. One of the pioneered methods from this group is the face recognition system, coined tensorfaces, proposed by Vasilescu and Terzopoulos [START_REF] Vasilescu | Multilinear analysis of image ensembles: TensorFaces[END_REF]. In their approach tensors are proposed to cope with multiple factors of face patterns, such as different poses, views, illuminations, etc. Another tensor based method for handwritten digits recognition was proposed by Savas et al. [START_REF] Savas | Handwritten digit classification using higher order singular value decomposition[END_REF] [START_REF] Lecun | Gradient-Based Learning Applied to Document Recognition[END_REF]. Their method assumes tensor decomposition which allows representation of a tensor as a product of its core tensor and a set of unitary mode matrices. This decomposition is called Higher-Order Singular Value Decomposition (HOSVD) [START_REF] Cichocki | Nonnegative Matrix and Tensor Factorization[END_REF][14] [START_REF] Kolda | Tensor Decompositions and Applications[END_REF]. A similar approach was undertaken by Cyganek in the system for road signs recognition [START_REF] Cyganek | An Analysis of the Road Signs Classification Based on the Higher-Order Singular Value Decomposition of the Deformable Pattern Tensors, Advanced Concepts for Intelligent Vision Systems Acivs[END_REF]. In this case, the input pattern tensor is built from artificially generated deformed versions of the prototype road sign exemplars. All aforementioned systems, which are based on HOSVD, show very high accuracy and high speed of response. However, computation of the HOSVD from large size tensors is computationally demanding since the algorithm requires successive computation of the SVD decompositions of matrices obtained from tensor flattening in different modes [START_REF] Lathauwer | Signal Processing Based on Multilinear Algebra[END_REF]. In practice, these matrices can be very large since they correspond to the products of all dimensions of the input tensor. In many applications this can be very problematic. To overcome this problem an ensemble with smaller size pattern tensor was proposed by Cyganek [START_REF] Cyganek | Ensemble of Tensor Classifiers Based on the Higher-Order Singular Value Decomposition[END_REF]. In the proposed methods tensors are of much smaller size than in a case of a single classifier due to the bagging process. However, despite the computational advantages, the proposed ensemble based method shows better accuracy when compared to a single classifier.

In this paper two modifications to the previously presented method are proposed. The first one is embedding of the Extended Euclidean Distance metric, recently introduced by Wang et al. [START_REF] Wang | On the Euclidean Distances of Images[END_REF]. This allows robustness to small geometrical perturbations of the input patterns since the new metric accounts for the spatial relationship of pixels in the input images. The second improvement consists in application of the weighted majority voting for combination of the responses of the classifiers in the ensemble. The experimental results show that in many cases the proposed improvements allow an increase of the overall accuracy of classification.

The rest of the paper is organized as follows. In section 2 properties of the Euclidean Image Distance are presented. In Section 3 the architecture of the proposed ensemble of the HOSVD multi-classifiers is discussed. Pattern recognition by the ensemble of the tensor classifiers is discussed in section 4. Experimental results are presented in section 5. The paper ends with conclusions in section 6.

Embedding Euclidean Image Distance

Images are 2D structures in which a scalar, vector (color) or multi-dimensional (MRI) value of a pixel is as important as its position within image coordinate space. However, the second aspect is not easy to be accounted for due to geometrical transformation of images of observed objects. On the other hand, image recognition heavily relies on comparison of images for which the Euclidean metric is the most frequently used one, mostly due to its popularity and simplicity in computations. However, Wang et al. proposed a better metric than Euclidean which takes into account also spatial relationship among pixels [START_REF] Wang | On the Euclidean Distances of Images[END_REF]. The proposed metric, called IMage Euclidean Distance (IMED), shows many useful properties, among which the most important is its insensitivity to small geometrical deformations of compared images.

More specifically, instead of the Euclidean metric between the two images X and Y of dimensions M×N each, given as follows

( ) ( ) ( ) ( ) = = - = - - ∑ 2 1 , MN T k k E k D x y X Y x y x y , (1)
Wang et al. propose to use the following extended version

( ) ( )( ) ( ) ( ) = = - - = - - ∑ , 1 , MN T k k l l G k l k l D g x y x y X Y x y G x y , (2)
where x and y are column vectors formed by the column-or row-wise vectorization of the images X and Y, respectively, and g kl are elements of the symmetric nonnegative matrix G of dimensions MN×MN, which defines the metric properties of the image space. Thanks to the above formulation, information on spatial position of pixels can be embedded into the distance measure, through the coefficients g kl . In other words, the closer the pixels are, the higher value of g kl should be, reaching its maximum for k=l. The distance between pixel positions (not values) is defined on an integer image lattice simply as a function of the 'pure' Euclidean distance between the points, as follows

( ) σ πσ - - = - = 2 2 2 2 1 2 k l kl k l g f e P P P P , (3)
where

⎡ ⎤ = ⎣ ⎦ 1 2
,

T i i i
p p P denotes position of the i-th pixel in the image, while σ is a width parameter, usually set to 1 [START_REF] Wang | On the Euclidean Distances of Images[END_REF]. Finally, incorporating (3) into (2) the IMED distance among image X and Y is obtained, as follows

( ) ( ) ( ) ( )( ) σ πσ - + - - = = - - ∑ 2 2 1 1 2 2 2 2 2 , 1 1 , 2 k l k l p p p p MN k k l l IMED k l D e x y x y X Y . ( 4
)
The D IMED image metric given in (4) can be used for a direct comparison of images, such as in the case of the k-nearest neighbor method, etc. It can be also incorporated into other classification algorithms, such as the discussed HOSVD. This can be achieved substituting D IMED into all places in which the D E was used.

However, for large databases of images direct computation of (4) can be expensive. An algorithm to overcome this problem was proposed by Sun et al. after observing that computation of D IMED can be equivalently stated as a transform domain smoothing [START_REF] Sun | A Fast Algorithm for Image Euclidean Distance[END_REF]. They developed the Convolution Standardized Transform (CST) which approximates well the D IMED . For this purpose the following separable filter was used

= ⊗ = T T H h h hh , (5)
where ⊗ denotes the Kronecker product of two vectors h and h T with the following components

⎡ ⎤ = ⎣ ⎦ 0.0053 0.2171 0.5519 0.2171 0.0053 T h . ( 6
)
The filter h given by ( 6) was also used in our computations since it offers much faster computations than direct application of (4).

Fig. 1 shows examples of application of the Standardizing Transformation for selected pictograms of the road signs in implementation with the filter h in [START_REF] Fukunaga | Introduction to Statistical Pattern Recognition[END_REF]. It is visible that ST operates as a low-pass filter (lower row). This way transformed patterns are fed to the classifier system. Finally, it should be noticed that the IMED transformation should not be confused with the Mahalanobis distance or the whitening transformation [START_REF] Fukunaga | Introduction to Statistical Pattern Recognition[END_REF] [START_REF] Duda | Pattern Classification[END_REF]. Specifically, in equation (2) we do not assume computation of any data distribution nor probabilistic spaces. In other words, the main difference lies in definition of the matrix G in (2) which elements, given by (3), convey information on mutual positions of the points. In contrast, for the Mahalanobis distance G would be an inverse of the covariance matrix which elements are computed directly from the values of x and y disregarding their placement in the images.

Architecture of the Ensemble of HOSVD Multi-Classifiers

Multidimensional data are handled efficiently with help of the tensor based methods since each degree of freedom can be represented with a separate index of a tensor [START_REF] Cichocki | Nonnegative Matrix and Tensor Factorization[END_REF] [START_REF] Cyganek | An Introduction to 3D Computer Vision Techniques and Algorithms[END_REF]. Following this idea, multidimensional training patterns can be efficiently represented by a prototype tensor [START_REF] Cyganek | An Analysis of the Road Signs Classification Based on the Higher-Order Singular Value Decomposition of the Deformable Pattern Tensors, Advanced Concepts for Intelligent Vision Systems Acivs[END_REF]. For the purpose of pattern recognition the prototype patterns tensor can be further decomposed into the orthogonal components which span a prototype tensor space. For the decomposition the Higher-Order Singular Value Decomposition can be used [1][13][11]. This way obtained orthogonal bases are then used for pattern recognition in a similar way to the standard PCA based classifiers [START_REF] Duda | Pattern Classification[END_REF][20] [START_REF] Turk | Eigenfaces for recognition[END_REF].

The HOSVD method allows any P-dimensional tensor

1 2 m n P N N N N N × × × × ∈ ℜ … … …
T to be equivalently represented in the following form [START_REF] Lathauwer | Signal Processing Based on Multilinear Algebra[END_REF][14]

1 1 2 2 P P = × × × S S S … T Z , (7)
where S k are N k ×N k unitary mode matrices,

1 2 m n P N N N N N × × × × ∈ ℜ … … …
Z is a core tensor. Z fulfills the following properties [START_REF] Lathauwer | Signal Processing Based on Multilinear Algebra[END_REF][14]:

1. (Orthogonality) Two subtensors

k k n a n b = = ⋅ = Z Z . ( 8
)
2. (Energy) All subtensors of Z for all k can be ordered according to their Frobenius norms, as follows

1 2 0 k k k P n n nN = = = ≥ ≥ ≥ ≥ … Z Z Z , (9)
The a-mode singular value of T is defined as follows In the next step, accuracies of each of the classifiers in the ensemble are assessed using the whole training dataset. These are then used to compute the weights of the classifiers in the ensemble. In the run-time their outputs are combined with the weighted majority voting scheme, as will be described in the next section.

Pattern Recognition by the Ensemble of Tensor Classifiers

It can be observed that thanks to the commutative properties of the k-mode tensor multiplication [START_REF] Savas | Handwritten digit classification using higher order singular value decomposition[END_REF], the following sum can be constructed for each mode matrix S i in ( 7)

1 P N h h P P h = = × ∑ s T T . ( 11
)
In the above the tensors

1 1 2 2 1 1 h P P - - = × × × S S S … T Z ( 12
)
form the basis tensors, whereas s h P denote columns of the unitary matrix S P . Because each T h is of dimension P-1 then × P in ( 11) is an outer product, i.e. a product of two tensors of dimensions P-1 and 1. Moreover, due to the orthogonality properties (8) of the core tensor Z in [START_REF] Kuncheva | Combining Pattern Classifiers[END_REF], T h are also orthogonal. Hence, they can constitute a basis which spans a subspace. This property is used to construct a HOSVD based classifier.

In the tensor space spanned by T h , pattern recognition can be stated as a measuring a distance of a given test pattern P x to its projections into each of the spaces spanned by the set of the bases T h in [START_REF] Kuncheva | Combining Pattern Classifiers[END_REF]. This can be written as the following minimization problem [START_REF] Savas | Handwritten digit classification using higher order singular value decomposition[END_REF]

] = -∑ 2 , 1 min i h i H i i x h h i c h Q c P T , ( 13
)
where the scalars c i h denote unknown coordinates of P x in the space spanned by T h i , H≤N P denotes a number of chosen dominating components.

To solve [START_REF] Lathauwer | Signal Processing Based on Multilinear Algebra[END_REF] the squared norm Q of ( 13) is created for a chosen index i. ˆ1 ,

Assuming further that

H i i h x h P T . ( 14
)
Thus, to minimize (13) the following value needs to be maximized

ρ = = ∑ 2 1
ˆ,

H i i h x h P T , (15)
Thanks to the above, the HOSVD classifier returns a class i for which its ρ i from [START_REF] Lecun | Gradient-Based Learning Applied to Document Recognition[END_REF] is the largest.

HOSVD0 p00 p00 p00 … HOSVD1 p10 p11 p12 … HOSVD2 p20 p21 p22 … … … … … …
In this work also different fusion methods were tested. Especially, the majority voting scheme was substituted for the weighted majority vote [START_REF] Kittler | On Combining Classifiers[END_REF] [START_REF] Polikar | Ensemble Based Systems in Decision Making[END_REF]. As alluded to previously, we proposed to use bagging to train the HOSVD classifiers from the ensemble which allows efficient memory usage. However, the partitions used for bagging contain less exemplars than all available for each prototype pattern. Therefore we further propose to use the whole training dataset to test each classifier trained with only fraction of that dataset for recognition of each pattern. This way we can assign some weight accuracies p kl for each classifier k and for each trained class l. These are defined as follows

= + l TP kl l l TP FP N p N N , (16)
where l TP N denotes a number of true positive responses and l FP N false positives, respectively. Table 1 visualizes this process for ten training patterns, such as digits. Further, it is assumed that the classifiers are independent and each is endowed with its individual accuracies p kl . If their outputs are combined with the weighted majority voting scheme, then accuracy of the ensemble is maximized by assigning the following weights [START_REF] Kuncheva | Combining Pattern Classifiers[END_REF]

[9] = - log 1 kl kl kl p b p , (17)
where p kl are given by ( 16). On the other hand, for a test pattern P x each of the HOSVD classifiers in the ensemble responds with its class and assigned vote strength, as follows

ρ ⎧ ⎪ = ⎨ ⎪ ⎩ ˆ, 0, kl k kl if HOSVD labels class l d otherwise , (18)
where ρ ˆkl denotes a maximal value of ρ ˆi in [START_REF] Lecun | Gradient-Based Learning Applied to Document Recognition[END_REF] and for the l-th classifier in ensemble and for pattern class k=i. Finally, the following discriminating function is computed

( ) ( ) = = + ∑ 1 ˆlog L k x k k l k l l g P db P , (19)
where P k denotes the prior probability for the k-th class. However, the latter is usually unknown, so in the rest of experiments the first term in [START_REF] Szeliski | Computer Vision. Algorithms and Applications[END_REF] was set to 0.

Experimental Results

The presented method was implemented in C++ using the HIL library [START_REF] Cyganek | An Introduction to 3D Computer Vision Techniques and Algorithms[END_REF].

Experiments were run on the computer with 8 GB RAM and Pentium® Quad Core Q 820 (clock 1.73 GHz).

For the experiments the USPS dataset was used [START_REF] Hull | A database for handwritten text recognition research[END_REF][24]. The same set was also used by Savas et al. [START_REF] Savas | Handwritten digit classification using higher order singular value decomposition[END_REF], as well as in the paper [START_REF] Cyganek | Ensemble of Tensor Classifiers Based on the Higher-Order Singular Value Decomposition[END_REF]. This dataset contains selected and preprocessed scans of the handwritten digits from envelopes of the U.S. Postal Service. Fig. 3 depicts some digits from the training and from the testing sets, respectively. The dataset is relatively difficult for machine classification since the reported human error is 2.5%. Therefore it has been used for comparison of different classifiers [START_REF] Lecun | Gradient-Based Learning Applied to Document Recognition[END_REF] [START_REF] Savas | Handwritten digit classification using higher order singular value decomposition[END_REF]. Originally the test and train patterns from the ZIP database come as the 16×16 gray level images. modifications proposed in this paper (the lower, the better). Blue bars for reference settings from [START_REF] Cyganek | Ensemble of Tensor Classifiers Based on the Higher-Order Singular Value Decomposition[END_REF]. Red bars relate to the method proposed in this paper, i.e. with the IMED and weighted majority voting.

Each experimental setup was run number of times and an average answer is reported. In all cases the Gaussian noise was added to the input image at level of 10%, in accordance with the procedure described in [START_REF] Cyganek | An Introduction to 3D Computer Vision Techniques and Algorithms[END_REF]. An analysis of different combinations of data partitions, number of classifiers in ensemble, number of dominating components, as well as input image size is discussed in our previous work [START_REF] Cyganek | Ensemble of Tensor Classifiers Based on the Higher-Order Singular Value Decomposition[END_REF]. In this paper we choose the best settings described in experimental results of the mentioned paper [START_REF] Cyganek | Ensemble of Tensor Classifiers Based on the Higher-Order Singular Value Decomposition[END_REF] and tested influence of the new IMED based preprocessing method, as well as new output fusion method. Results for three different settings are shown in the bar graph in Fig. 4.

For the experimental setups in this work were chosen three best setups from our previous work [START_REF] Cyganek | Ensemble of Tensor Classifiers Based on the Higher-Order Singular Value Decomposition[END_REF]. These are summarized in Table 2. The first experiments were run with configurations of the ensembles from Table 2. Then new propositions of IMED and weighted majority voting were introduced and run again. Each experiment setup was run 25 times and average parameters are reported. In all cases the proposed modifications allowed better results of 0.2-1% with negligible time penalty due to separability of the IMED filter.

Conclusions

In this paper an extended version of our previous work on image classification with the ensemble of tensor based classifiers is presented [START_REF] Cyganek | Ensemble of Tensor Classifiers Based on the Higher-Order Singular Value Decomposition[END_REF]. In this method, thanks to the construction of the ensemble of cooperating classifiers, tensors are of much smaller size than in a case of a single classifier. Each classifier in this ensemble is trained with data partition obtained from bagging. Such approach allows computations with much smaller memory requirements. However, the method shows also better accuracy when compared to a single classifier. In this paper two modifications to this formulation were discussed. The first is to apply input pattern preprocessing with embedding of the Extended Euclidean Distance metric. This allows robustness to small geometrical perturbations thanks to the metric which accounts for the spatial relationship of pixels in the input images. The second improvement consists in application of the weighted majority voting for combination of the responses of the classifiers in the ensemble. The experimental results show that in many cases the proposed improvements allow an increase of the overall accuracy of classification in order of 0.2-1%. The method is highly universal and can be used with other types of patterns.

Fig. 1 .

 1 Fig. 1. Visualization of Standardizing Transform applied to the road sign pictograms. Original pictograms (upper row). After transformation (lower row).

Z

  for all possible values of k for which a≠b it holds that 0

  computation of the HOSVD is based on successive computations of the SVD decomposition of the matrices composed of the flattened version of the tensor T . The algorithm requires a number of the SVD computations which is equal to the valence of that tensor. The detailed algorithm can be referred to in the literature [1][13][14].

Fig. 2 .

 2 Fig. 2. Architecture of the ensemble with the HOSVD classifiers. Data preprocessed with the Image Euclidean Transformation. Bagging method used for training. Outputs are combined with the weighted majority voting.

Fig. 2

 2 Fig. 2 shows architecture of the proposed ensemble of the HOSVD classifiers. All training and testing data are preprocessed with the Image Euclidean Transformation

T

  and x P are normalized the following is obtained (the hat mark indicates tensor normalization)

Fig. 3 .Fig. 4 .

 34 Fig. 3. Two data sets from the ZIP database. Training set (a), and testing set (b) The database is divided into the training and testing partitions, counting 7291 and 2007 exemplars, respectively.

Table 1 .

 1 Structure of a matrix of partial accuracies for each classifier and each training prototype pattern

	Digit	0	1	2	3	4	5	6	7	8	9

Table 2 .

 2 Three best experimental setups of the ensemble from[START_REF] Cyganek | Ensemble of Tensor Classifiers Based on the Higher-Order Singular Value Decomposition[END_REF]

	No.	Param.	Number of experts	Data in samples	Important components	Image resolution	Noise [%]
	1		11	64	16	16x16	10
	2		15	192	16	32x32	10
	3		33	64	16	16x16

Acknowledgement

The work was supported in the years 2011-2012 from the funds of the Polish National Science Centre NCN, contract no. DEC-2011/01/B/ST6/01994.


Václav Snášel
email: vaclav.snasel@vsb.cz

Pavel Krömer
email: pavel.kromer@vsb.cz

Jan Platoš
email: jan.platos@vsb.cz

Miloš Kudělka
email: milos.kudelka@vsb.cz

Zdeněk Horák
email: zdenek.horak@vsb.cz

Floriana Esposito

Attilio Giordana

Lorenza Saitta

Ugo Galassi

Alessandro Serra

Giuseppe Lo Bello

M E Campidoglio -

Filippo Neri

Davide De

Cristina Baroglio

Filippo Brancadori

M Paderni -

R Finelli

Marco Botta

Davide Roverso

Pietro Laface

Renato De Mori

On Spectral Partitioning of Co-authorship Networks

Keywords: spectral partitioning, algebraic connectivity, co-authorship, DBLP

Spectral partitioning is a well known method in the area of graph and matrix analysis. Several approaches based on spectral partitioning and spectral clustering were used to detect structures in real world networks and databases. In this paper, we explore two community detection approaches based on the spectral partitioning to analyze a coauthorship network. The partitioning exploits the concepts of algebraic connectivity and characteristic valuation to form components useful for the analysis of relations and communities in real world social networks.

Introduction

Spectral clustering (or spectral partitioning) is a useful method for partitioning and clustering of graphs and networks with solid mathematical background and clear interpretation. The ubiquity of social and communication networks in today's information society hand in hand with the increasing power of computers makes the usage of algebraic techniques such as spectral clustering very practical. In this work, we use the spectral partitioning to analyze selected parts of the DBLP1 , a large database of computer science publications. The DBLP can be seen as a vast, dynamic and constantly updated social network that captures several years of author co-operations in the form of joint publications. It is very interesting for social network (SN) researcher because the authors can be easily grouped based on their affiliations, areas of interest, and advisor-advisee relationship. Moreover, we can trace in the DBLP the development of each author's activities, types of activities, areas of interest and so on.

In this paper, we present two spectral partitioning based algorithms to iteratively detect communities in the DBLP (and social networks in general).

The basics of the spectral clustering (SC) were introduced in 1975 by M. Fiedler [START_REF] Fiedler | A property of eigenvectors of nonnegative symmetric matrices and its application to graph theory[END_REF]. Fiedler's work defined spectral clustering for both, unweighted and weighted graphs. The following definitions apply to weighted graphs because the edges in a co-authorship network intuitivelly have different weights. An edge between two authors that have published one joint paper has different quality (i.e. weight) than an edge between two authors that have published a large number of joint papers through the years. The frequency, regularity, and age of such co-operations can be a hint for an edge weighting scheme.

Definition 1 (Generalized Laplacian of weighted graph G) For a graph G = (V, E), the generalized Laplacian is the matrix of the quadratic form

(A C (G)x, x) = (i,k)∈E c ik (x i -x k ) 2 (1)
The A C (G) can be easily computed:

a ik = 0 if i = k and (i, k) ∈ E -c ik if i = k and (i, k) ∈ E (2)
a ii = - k =i c ik i, k ∈ N (3)
Definition 2 (Algebraic connectivity of weighted graph G) The algebraic connectivity of graph G denoted a C (G) is the second smallest (first non-zero)

eigenvalue of A C (G). Let 0 = λ 1 ≤ λ 2 ≤ . . . ≤ λ n be the eignevalues of A C (G). Then a C (G) = λ 2 .
The algebraic connectivity a C (G) is also known as the Fiedler value [START_REF] Spielman | Spectral partitioning works: Planar graphs and finite element meshes[END_REF].

Definition 3 (Characteristic valuation of G)

The characteristic valuation of G (also known as the Fiedler vector of G) denoted a(G) = (a 1 , . . . , a n ) is defined by the values of the eigenvector corresponding to a C (G).

The characteristic valuation assigns a non-zero (positive or negative) value to each vertex in the graph in a natural way. There is a number of interesting properties of a C (G) and a, for example [START_REF] Fiedler | A property of eigenvectors of nonnegative symmetric matrices and its application to graph theory[END_REF][START_REF] Grady | Discrete Calculus -Applied Analysis on Graphs for Computational Science[END_REF][START_REF] Spielman | Spectral partitioning works: Planar graphs and finite element meshes[END_REF]:

-a C (G) is positive iff G is connected. -if a C (G)
is small, then a graph cut according to the values of vertices in a(G) will generate a cut with good ratio of cut edges to separated vertices. a(G) represents an ordering (Fiedler ordering) which can be used for spectral partitioning of connected graphs (for the rationale see theorem 1).

Theorem 1 For a finite connected graph G with n vertices that has a positive weight c ik assigned to each edge (i, k), characteristic valuation a(G), and any r ≥ 0 let

M (r) = {i ∈ N |y i + r ≥ 0} (4)
The subgraph G(r) induced by G on M (r) is connected.

Via theorem 1 can be defined iterative (stepwise) partitioning of connected graph G into connected subgraph G(r) and general subgraph G \ G(r). Via theorem 1 can be also defined iterative elimination of vertices with lowest significance to the graph so that the remainder of the graph is connected. The proof of theorem 1 can be found in [START_REF] Fiedler | A property of eigenvectors of nonnegative symmetric matrices and its application to graph theory[END_REF].

Graph partitioning

A graph G = (V, E) can be partitioned into two disjoint sets A, B such that A ∪ B = V and A ∩ B = ∅. The cut value, which describes the dissimilarity between the two partitions, can be defined as the sum of weights of the edges removed by the cut [START_REF] Shi | Normalized cuts and image segmentation[END_REF]:

cut(A, B) = i∈A,j∈B c ij (5)
It can be shown that the Fiedler vector represents solution for finding partitions A and B such that the following cost function (the average cut) is minimized [START_REF] Shen | Graph-theory based parcellation of functional subunits in the brain from resting-state fmri data[END_REF][START_REF] Shi | Normalized cuts and image segmentation[END_REF]:

Acut(A, B) = cut(A, B) |A| + cut(B, A) |B| (6)
The average cut is a measure with known imperfections [START_REF] Shi | Normalized cuts and image segmentation[END_REF]. However, its usage is simple and its computation is fast.

Related work

As the need for efficient analysis of graph-like structures including social networks is growing, there was much attention given to spectral partitioning and spectral clustering of graphs. In this section, we provide brief state of the art of graph partitioning methods based on spectral clustering. The use of spectral partitioning for graph analysis was advocated by Spielman and Teng [START_REF] Spielman | Spectral partitioning works: Planar graphs and finite element meshes[END_REF]. They have shown that spectral partitioning works well for bounded-degree planar graphs and well-shaped d-dimensional meshes. Today, methods based on spectral clustering are being used to analyze the structure of a number of networks.

An influential study on spectral clustering and its application to image segmentation was published in 2000 by Shi and Malik [START_REF] Shi | Normalized cuts and image segmentation[END_REF]. The authors approached the graph partitioning task from the graph cuts point of view. They described the graph cut defined by the Fiedler vector and called it average cut. The average cut was shown to be good at finding graph splits whereas the newly defined normalized cut was designed to compute the cut costs as a ratio of cut edge weights to all edge weights in the segments. The normalized cut was shown to be useful when seeking partitions that are both, balanced and tight. On the other hand, a study by Sarkar and Soundararajan showed that the increased computational cost of the normalized cut does not result in statistically better partitions [START_REF] Sarkar | Supervised learning of large perceptual organization: graph spectral partitioning and learning automata[END_REF].

Ding et al. [START_REF] Ding | A min-max cut algorithm for graph partitioning and data clustering[END_REF] have proposed in 2001 another graph cut algorithm, the minmax cut, and showed its usefulness for partitioning real world graphs into balanced parts. Bach and Jordan [START_REF] Bach | Learning spectral clustering[END_REF] proposed an algorithm based on a new cost function evaluating the error between given partition and a minimum normalized graph cut. The partitions can be learned from given similarity matrix and vice-versa -the similarity matrix can be learned from given clusters. Similarity of nodes i and j in this context means large weight of the edge (i, j), i.e. large c ij . The method leads to clusters with large in-cluster similarity and small inter-cluster similarity of nodes.

The algebraic connectivity has been used to define a new method for construction of well-connected graphs by Gosh and Boyd in 2006 [START_REF] Ghosh | Growing well-connected graphs[END_REF]. The algorithm uses the properties of algebraic connectivity and defines an edge perturbation heuristic based on the Fiedler vector to choose from the set of candidate edges such edges that would improve the value of a C (G).

The work of Ruan and Zhang [START_REF] Ruan | An efficient spectral algorithm for network community discovery and its applications to biological and social networks[END_REF] presents an application of spectral partitioning in the area of social networks. The authors developed an efficient and scalable algorithm Kcut to partition the network to k components so that the modularity Q of community structures is maximized. For more details on Q see [START_REF] Ruan | An efficient spectral algorithm for network community discovery and its applications to biological and social networks[END_REF]. The usefulness and effectiveness of Kcut was demonstrated on several artificial and real world networks.

Mishra et al. [START_REF] Mishra | Clustering social networks[END_REF] have used spectral clustering for social network analysis in 2007. They aimed at finding good cuts on the basis of conductance, i.e. the ratio of edges crossing the cut to the minimum volume of both partitions. Volume in this context means the number of edges incident with vertices in the sub-graph. Moreover, the proposed algorithm was able to find overlapping clusters with maximum internal density and external sparsity of the edges. Kurucz et al. [START_REF] Kurucz | Spectral clustering in telephone call graphs[END_REF][START_REF] Kurucz | Advances in web mining and web usage analysis[END_REF] have applied spectral clustering to telephone call graphs and to social networks in general. In their studies, the authors discussed various types of Laplacians, edge weighting strategies, component size balancing heuristics, and the number of eigenvectors to be utilized. The work proposed a k-way hierarchical spectral clustering algorithm with heuristic to balance clusters and showed its superiority over the Divide-and-Merge clustering algorithm.

In 2008, Leskovec et al. [START_REF] Leskovec | Statistical properties of community structure in large social and information networks[END_REF] investigated the statistical properties of communities in social and information networks. They used the network community profile plot to define communities according to the conductance measure. Their work demonstrated that the largest communities in many real world data sets blend with the rest of the graph with increasing size, i.e. their conductance score is decreasing.

Xu et al. [START_REF] Xu | Revealing social networks of spammers through spectral clustering[END_REF] have analyzed social networks of spammers by spectral clustering. They have used the normalized cut diassociation measure that is known to minimize the normalized cut between clusters and simultaneously maximize the normalized association within clusters.

A recent work on generalized spectral clustering based on the graph p-Laplacian is due to Bühler and Hein [START_REF] Bühler | Spectral clustering based on the graph p-laplacian[END_REF]. It was shown that for p → 1 the cut defined by Fiedler vector converges to the Cheeger cut. The p-Spectral Clustering using the p-Laplacian, a nonlinear generalization of the graph Laplacian, was in this paper evaluated on several data sets.

An overview of spectral partitioning with different Laplacians was given by Luxburg in [START_REF] Luxburg | A tutorial on spectral clustering[END_REF]. The study contained a detailed description of the algorithm, properties of different Laplacians and a discussion on suitability of selected Laplacians for given task.

In general, many variants of the basic spectral clustering algorithm were used to partition graphs and detect network structure in multiple application areas with good results. Real world networks and social networks constituted by the natural phenomena of communication, interaction, and cooperation are especially interesting application field for the spectral partitioning.

Spectral partitioning of co-author communities in the DBLP

We have defined two iterative partitioning algorithms based on spectral clustering and algebraic connectivity to find co-author communities in the graph. In the algorithm 1 (simple iterative spectral partitioning, SimpleISP) was the initial connected graph divided into two subgraphs, each containing vertices with positive valuation (and incident edges) and vertices with negative valuation (and incident edges) respectively. For the next iteration was used as an input the subgraph that contained the author vertex. If the author vertex belonged to the negative subgraph (that was not guaranteed to be connected), all vertices that were not connected to the author vertex were removed. The partitioning ended when the subgraph contained only single vertex (author vertex ). This variant of the algorithm creates in every iteration a smaller (narrower) community centered around the author.

In the algorithm 2 (iterative spectral patritioning, ISP), the graph for next iteration was created differently. In each iteration, we removed all vertices that had lower characteristic valuation than the author vertex. It is guaranteed that the resulting subgraph is connected. The algorithm ended when the author vertex had the lowest valuation among all vertices in the graph (i.e. it was not possible to remove loosely connected vertices). This variant of the algorithm centers on the community to which the author belongs rather than on the author herself.

Algorithm 1 Simple iterative spectral partitioning (SimpleISP)

1: Find a connected subgraph S containing the vertex of selected author (author vertex ), vertices of all his or her co-authors, vertices of all their co-authors, and edges among them.

2: while |S| > 1 do 3:

Compute a(S)

4:

Cut S according to a(S)

5:

Let S + contain all vertices and incident edges for which the value of a(S)i ≥ 0 and S - contain all vertices and incident edges for which a(S)i < 0.

6:

Remove all edges between vertices in S + and S -.

7:

if author vertex ∈ S + then 8:

S = S + 9:

12:

Remove from S all vertices that are not connected to author vertex 13: end while Algorithm 2 Iterative spectral partitioning (ISP)

1: Find a connected subgraph S containing the vertex of selected author (author vertex ), vertices of all his or her co-authors, vertices of all their co-authors, and edges among them.

2: repeat 3:

Compute a(S)

4:

Get the valuation of author vertex a AV = a(S) author vertex

5:

Remove from S all vertices with valuation lower than a AV . The rest is connected.

6: until min a(S) < a AV

Experiments

To observe the communities generated by proposed algorithms, we have conducted a series of experiments with the DBLP data. We have downloaded the DBLP dataset from April 2010 in XML and preprocessed it for further usage. We have selected all conferences held by IEEE, ACM or Springer, which gave us 9,768 conferences. For every conference we identified the month and year of the conference. In the next step we extracted all authors having at least one published paper in the mentioned conferences (as authors or co-authors). This gave us 443,838 authors. Using the information about authors and their papers we were able to create a set of cooperations between these authors consisting of 2,054,403 items. Finally, the cooperations were represented as a graph. A vertex in the graph represented one author and an edge represented a co-operation between the authors (joint publication). The edges were weighted according to the number of joint publications between the two authors, i.e. if two authors published one joint work, the weight of the edge between their vertices was 1. If they co-operated on n papers, the weight of the edge between their vertices was n. We note that this weighting scheme is quite naïve and much more sophisticated approaches can be used, but such a research is out of the scope of this paper.

We have selected two authors and investigated spectral partitions of the connected graph consisting of their co-authors and their co-authors' co-authors.

We investigated only two levels of co-authors to obtain components that could be manually inspected. Floriana Esposito and Philip S. Yu were investigated in a recent work on co-authorship network analysis [START_REF] Kudělka | Social and swarm aspects of co-authorship network[END_REF]. Floriana Esposito is an author who has been active since 1990 and who has a lot of strong ties whereas Philip S. Yu is an author with the greatest number of records in the data set and with a number of strong co-authors. We have applied both, simple iterative spectral partitioning and strict iterative spectral partitioning to the subgraphs around selected authors.

Results

The process of iterative spectral partitioning of subgraphs for Philip S. Yu and Floriana Esposito is captured in Fig. 1. The figures illustrate the sizes of components (communities) of both authors in each iteration of SimpleISP and ISP. Initial size of P. S. Yus component was 9607 and initial size of F. Espositos component was 1180, so for a better comparison, the relative component sizes are compared in Fig. 1(c) and Fig. 1(d). Figures Fig. 1(a) and Fig. 1(c) show the process of SimpleISP. We can see that both authors loose the majority of their collaborators in the second iteration. However, Floriana Espositos network keeps larger fraction of the original nodes during the whole process and it be- comes larger than Philip S. Yus network after 12th iteration. The SimpleISP ended for Floriana Esposito after 20 iterations and for Philip S. Yu after 16 iterations.

The ISP process is shown in Fig. 1(b) and Fig. 1(d). In this case, the most significant reduction of the communities was done in the first iteration. Floriana Espositos network lost 932 out of 1180 nodes and Philip S. Yus network reduced from 9607 to 410 nodes. Again, the relative component size of Floriana Esposito was greater than the relative component size of Philip S. Yu during the whole ISP and it becomes larger than P. S. Yus community after second iteration. The ISP ended for F. Esposito after 7 iterations and the final network conatined 20 nodes. In contrast, the ISP for P. S. Yu ended after 8 iterations and the final network contained only one node -the author node.

Examples of the partitions in selected iterations of the SimpleISP and ISP for P. S. Yu and F. Esposito are shown in Fig. 2, Fig. 3, Fig. 4, and Fig. 5 respectively. Blue and red vertices and edges represent the components and dotted edges represent the cut. The number on each vertex corresponds to characteristic valuation of the vertex and the number on each edge represents the weight of the edge, i.e. the multiplicity of author co-operation in this experiment. We note that larger graphs are shown to illustrate the structure of the community and cut rather than to provide the names of the co-authors which is printed using very small font.

Conclusions and future work

In this paper we present two algorithms for iterative spectral partitioning of social networks. The goal of the algorithms is to find meaningful communities in networked data. We demonstrate the application of the algorithms on a coauthorship network, namely the DBLP, in which we sought for communities of selected authors. The first algorithm focused on a central node around which it iteratively created connected subgraphs, i.e. possible communities. It was searching for communities around an author. The second algorithm, following more closely the idea of algebraic connectivity and spectral clustering, focused on a community rather than on the author. It was highlighting the community to which the author belonged. We have selected two authors with different statistical properties and searched for their communities using both approaches. The results of the experiment show that both, the partitioning process and generated partitions, were quite different for the two authors, no matter which algorithm was used. An author with strong ties to other authors retained connection to a large number of co-author nodes during most of the partitioning process. On the other hand, a highly co-operative author lost the links to majority of his/her co-authors very early. The results support the intuition that the

  Component size (SimpleISP).

  Relative component size (ISP).

Fig. 1 :

 1 Fig. 1: Size of author components during the partitioning.

Fig. 3 :

 3 Fig. 3: Floriana Espositos network in selected iterations of SimpleISP.

http://www.informatik.uni-trier.de/ ley/db/

Acknowledgments. This work was supported by the European Regional Development Fund in the IT4Innovations Centre of Excellence project (CZ.1.05/ 1.1.00/02.0070) and by the Bio-Inspired Methods: research, development and knowledge transfer project, reg. no. CZ.1.07/2.3.00/20.0073 funded by Operational Programme Education for Competitiveness, co-financed by ESF and state budget of the Czech Republic.

partitioning of such a different authors will be different. We have also observed, that the author with strong relationships to others was placed to a community of twenty collaborators whereas the highly collaborative author ended alone. There are many directions in which this work can continue. First, the observations presented in this paper should be confirmed on a large number of authors. Second, the weighting scheme used in this study was rather simple -a different edge weighting schemes should be applied and their influence on the partitioning should be investigated. Third, in this work we have used the simple average cut in which we have split the network according to negative and positive values of vertex characteristic valuation. Many different cuts were proposed and their effect on co-authorship network partitioning should be investigated. Also the effect of different Laplacians should be investigated. Finally, the results of the spectral clustering of the co-authorship network should be compared to other non-spectral network and graph analytical methods.


Jiří Dvorský
email: jiri.dvorsky@upol.cz

Zbyněk Janoška
email: zbynek.janoska@cdv.cz

Lukáš Vojáček
email: lukas.vojacek@vsb.cz

P Systems for Traffic Flow Simulation

Keywords: membrane computing, P systems, traffic flow, XML

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

Introduction

Membrane computing represents new and rapidly growing branch of natural computing, which starts from observation that the processes taking place in a living cell can be understood as a computation. Membrane computing and its computational device -P system -were introduced by Pȃun [START_REF] Pȃun | Computing with Membranes[END_REF] and gained a lot of interest in last decade. P systems start from observation, that membrane plays a fundamental role in the functioning of a living cell. Membranes act as threedimensional compartments which delimit various regions of a living cell. They are essentially involved in a number of reactions taking place inside cell and moreover act as selective channels of communication between different compartments of a cell [START_REF] Bernardini | On P Systems as a Modelling Tool for Biological Systems[END_REF].

P systems take inspiration from cell on two levels -the structure and the functioning. Structure of cell is represented by its membranes and functioning is governed by biochemical reactions. Every P system therefore has three main elements: a membrane structure, where object evolve according to given evolution rules [START_REF] Pȃun | Introduction to Membrane Computing[END_REF]. Some authors add fourth basic element of membrane systems -communication [START_REF] Bernardini | On P Systems as a Modelling Tool for Biological Systems[END_REF][START_REF] Pȃun | The Power of Communication: P Systems with Symport/Antiport[END_REF]. Communication is always encoded in rules (they are called communication rules instead of evolution rules) and will be dealt with later in the text. From the point of view of transportation modelling, communication (e.g. topology) is essential feature. Fig. 1. Graphical representation of P system, [START_REF]P system -Wikipedia, the free encyklopedia[END_REF] Simple example of P system is depicted in Fig. 1. Membrane structure is hierarchically arranged set of membranes, contained in a distinguished outer membrane, called skin membrane. System is surrounded by the environment, which may collect objects leaving the system, or in some variants of P systems, the environments can actively support system with objects [START_REF] Bernardini | Dynamical aspects of P systems[END_REF][START_REF] Cavaliere | Evolution-Communication P Systems[END_REF]. Membranes delimit regions, with which they are in one-toone relation. Therefore the terms membrane and region are mostly interchangeable. Each membrane is identified by its label, which can be with membranes in one-to-many relation. The position of inner membranes does not matter; we assume, that in membrane there is no ordering, everything is close to everything else [START_REF] Pȃun | Introduction to Membrane Computing[END_REF].

Second basic element of P systems are objects. By objects in biological sense are meant chemicals, ions, molecules etc. Those substances are present in a cell in enormous amounts, but the ordering again does not matter. What matters is the concentration, the population, the number of copies of each molecule [START_REF] Pȃun | Introduction to Membrane Computing[END_REF]. Abstracting from biological reality, we represent each substance by a symbol from given alphabet and since the multiplicity matters, instead of objects we use multisets of objects. Common notation of multisets in P systems is following: if, for example, objects a, b, c are present in 7,2 and 5 copies, they will be represented by multiset a 7 b 2 c 5 .

In basic variant of P systems, multisets of objects are considered to be floating in inner regions of membrane systems. They evolve by the means of evolution rules, which are localised with the regions of the membrane structure. There are three main types of rules [START_REF] Pȃun | Introduction to Membrane Computing[END_REF]: (1) multiset-rewriting rules, (2) communication rules and (3) rules for handling membranes. In this section only first type of rules will be described.

Multiset-rewriting rules take form u → v, where u and v are multisets of objects. For example, rule ab → cd 2 says, that one copy of a and one copy of b are consumed and one copy of c and two copies of d are produced. A number of possible extensions of rules exists.

Two crucial features of P systems have to be mentioned at this point. As mentioned earlier, in membranes everything is close to everything else. Therefore, if one instance of an object can be processed by two or more rules, the rule to be applied is chosen non-deterministically. All rules have the same probability to be chosen. The rules also have to be used in maximally parallel manner.

More specifically, the objects are assigned to rules, non-deterministically choosing the objects and the rules, until no further assignment is possible. An evolution step in a given region of membrane system consists of finding the maximal applicable multiset of rules, removing from region all objects specified in the left hand of the chosen rules and producing the objects on the right hand side of the rules.

After giving short introduction to basic notions of P systems, let us continue with more detailed survey on traffic flow modelling and representation of transportation system by the means of P systems. In Sect. 2, a definition of basic model -transitional P system, is given. Sect. 3 presents possible extensions of such model, which are incorporated into a model, described in Sect. 5. Sect. 4 gives a short introduction into the traffic flow modelling. Finally, an XML scheme designed to store the configuration of developed system is presented in Sect. 6. We will conclude with some final remarks in Sect. 7.

Transition P system

P systems based on application of multiset-rewriting rules are called transition P system. Formally, transition P system is a construct of the form:

Π = (O, C, µ, w 1 , w 2 , . . . , w m , R 1 , R 2 , . . . , R m , i o ), (1)
where:

-O is the finite and non-empty alphabet of objects, -C ⊂ O is the set of catalysts, µ is a membrane structure, consisting of m membranes, labeled 1, 2, . . . , m; one says, that the membrane structure, and hence the system, is of degree m, w 1 , w 2 , . . . , w m are strings over O representing multisets of objects present in regions 1, 2, . . . , m of membrane structure, -R 1 , R 2 , . . . , R m is finite set of evolution rules associated with regions 1, 2, . . . , m of membrane structure, i o is either one of the labels 1, 2, . . . , m and then the respective region is the output region of the system, or it os 0 and then the result of the computation is collected in the environment of the system.

A sequence of transitions of P system constitutes a computation. A computation is successful if it halts, it reaches a configuration where no rule can be applied to the existing objects, and output region i o still exists [START_REF] Pȃun | Introduction to Membrane Computing[END_REF].

The rules are of form u → v, where u ∈ O and v ∈ (O × T ar), where T ar = {here, in, out}. Target indications T ar extend transition P system in following way: rule ab → c here d in e out consumes one instance of each a and b and produces one copy of c in current membrane, one copy of d in a child of current membrane and one copy of e in the parent of current membrane. If current membrane is skin membrane, object e is send to environment of the system. If current membrane does not have a child, rule can not be applied.

Another extension comes from the existence of catalysts. Catalysts are objects, which participate in a chemical reaction, but are not consumed or produced by it. They just enable the application of rule. Rule with catalysts takes following form: ac → bc, with object c being the catalyst.

Possible Extensions of P systems

In this section we will mention some elementary extensions of P system, which however constitute only a fracture of possibilities. We refer reader to The P systems Webpage [START_REF]The P Systems Web Page[END_REF] for complete list of publications and further information. Already in the text, three types of rules were mentioned. Evolution rules were briefly covered in previous sections.

Communication rules were introduced in [START_REF] Pȃun | The Power of Communication: P Systems with Symport/Antiport[END_REF]. Basic idea of communicating P systems is, that computation is achieved only by transporting object between membranes. Direct inspiration from biology are symport and antiport. When two chemical pass through membrane only together, in the same direction, the process is called symport. When the two chemicals pass only with help of each other, but in opposite directions, the process is called antiport [START_REF] Pȃun | The Power of Communication: P Systems with Symport/Antiport[END_REF]. Symport rules take a form (ab, in) or (ab, out), and antiport rules take a form (a, out; b, in), where a, b are object from alphabet of all possible objects. Meaning of rules is following: for symport rule (ab, in) or (ab, out), if objects a, b are present in current membrane, they are sent together into child (or parent, in second case) of the membrane. For antiport rule (a, out; b, in), if a is present in current membrane and b is present in the parent of a membrane, than a exits current membrane and b enters it. Universality of P systems with symport and antiport have been proven [START_REF] Pȃun | The Power of Communication: P Systems with Symport/Antiport[END_REF] and simplified version of communication, conditional uniport have been studied [START_REF] Verlan | Computational completeness of tissue p systems with conditional uniport[END_REF]. Comprehensive review of communication strategies in P systems can be found in [START_REF] Verlan | Generalized communicating P systems[END_REF].

Third type of rules are rules for handling membranes. Dissolution of membranes has already been mentioned, but other ways to obtain dynamical membrane structure, evolving during the course of computation, have been presented. Most simple of those is assigning electrical polarization +, -, 0 to each membrane. Polarization replaces target indicators in, out, here. Polarization of objects is introduced by the rules and polarized objects can enter only membranes with opposite polarization. For example, ab → c + d -means, that one instance of c enters inner membrane with negative polarization and one instance of d enters inner positive membrane. Rules can also be used to change the polarization of membranes during the computation.

Two issues seem essential, when P systems are used to simulate real-world phenomena rather than for computation. Non-determinism is first of the issues. Some chemical reactions are more likely to occur than others. First attempt to solve this is by assigning priorities to rules. Firstly, set of rules with the highest priority is chosen and according to the principle of maximal parallelism, all rules which can be applied, are applied. Then, the rules with second highest priority are selected and the procedure repeats, etc.

Second approach is to assign probabilities to all rules. Probability can be introduced to P system on different levels [START_REF] Obtulowicz | In search of) Probabilistic P systems[END_REF], but here we will mention only probability on the level of rule selection. Different approaches have been proposed [START_REF] Bernardini | Dynamical aspects of P systems[END_REF][START_REF] Cazzaniga | Stochastic Approaches in P Systems for Simulating Biological Systems[END_REF][START_REF] Pescini | Dynamical probabilistic P systems[END_REF]. Basic idea is to associate each rule with a constant k, so the rule takes a form u k -→ v, where u, v are multisets of objects and k can be interpreted either as a probability, or as a "stoichiometric coefficient", using which the true probability is calculated.

Last extension, which we will mention at this point, is representation of time of P systems. In real world, every biochemical reaction takes some time. Representation of time in P systems is similar to representation of probability. A constant t is assigned to each rule, so the rule takes the form u t -→ v, where u, v are multisets of objects and t is number of time units, which must pass to complete the application of the rule [START_REF] Cavaliere | Time-independent P systems[END_REF]. In the first time step, multiset u is consumed and removed from the current membrane. After t -1 more time steps, multiset v is introduced into the system. Time can also be introduced into P systems as a lifetime of objects or even membranes [START_REF] Aman | Adding Lifetime to Objects and Membranes in P Systems[END_REF].

For the sake of brevity, we will not discuss more extensions of P systems, although many possibilities were explored within this framework. P system based model of traffic flow will take advandage of three main features mentioned in this section -membrane polarization, probability and time-dependence.

Traffic flow simulation

Hoogendoorn and Bovy [START_REF] Hoogendoorn | State-of-the-art of Vehicular Traffic Flow Modelling[END_REF] distinguish three main levels of modelling of transportation systems:

-microscopic simulation model describes both the space-time behaviour of the systems' entities (i.e. vehicles and drivers) as well as their interactions at a high level of detail (individually). Sometimes submicroscopic models are mentioned separately. -mezoscopic model does not distinguish nor trace individual vehicles, but specifies the behaviour of individuals, for instance in probabilistic terms. To this end, traffic is represented by (small) groups of traffic entities, the activities and interactions of which are described at a low detail level. -macroscopic flow model describe traffic at a high level of aggregation as a flow without distinguishing its constituent parts. For instance, the traffic stream is represented in an aggregate manner using characteristics as flowrate, density, and velocity. Individual vehicle manoeuvres, such as a lanechange, are usually not explicitly represented [START_REF] Hoogendoorn | State-of-the-art of Vehicular Traffic Flow Modelling[END_REF].

While microscopic models describe accurately behaviour of small-scaled systems such as lanes or intersections, their use to simulate spatially extensive areas, such are cities or agglomerations, are limited. Macroscopic models, on the other hand, are capable of simulation on low-scale level, but lack detailed description of individuals. Mezoscopic models are a reasonable compromise between detail and robustness, which is desirable for certain applications, such are public transport modelling [START_REF] Peeta | Foundations of dynamic traffic assignment: The past, the present and the future[END_REF].

In this study we propose a simulation model for traffic flow, aiming to describe accurately vehicle behaviour at the intersections, but ignoring detailed description of cars behaviour at road segments between intersections. Such model should be individual based, representing vehicles as agents instead of populations, but at the same time should be robust enough to enable simulation of large-scaled road networks. Such model could be used for example for optimalization of traffic lights at signalized intersections.

P system based model of traffic flow

In this section we describe P system model of traffic flow. This model aims to capture detailed behaviour of cars at the intersection, while ignoring behaviour at road segments.

Definition of a system

Within P systems paradigm, road network can be described as a graph with membranes at the nodes. This representation is similar to the one of communication P systems [START_REF] Verlan | Generalized communicating P systems[END_REF]. There are three types of membranes, based on their function. White holes act as generators of cars. If the systems represents, say, a city, than white holes will be the roads, bringing cars to the city. Opposite reactions take place in black holes, where the cars leave the system. White and black holes does not necessarily have to be only at the periphery of the system -also residential, industrial or business areas within the city can act as white or black holes, and they can therefore be placed arbitrarily in the graph. The third type of nodes is an intersection. At intersection, cars are distributed to other parts of the network. Within this paradigm, cars are represented as objects. For the sake of brevity, we will consider only one type of object -an arbitrary vehicle. In praxis, however, different types of vehicles (cars, buses, motorcycles etc.) can be represented by different objects.

The behaviour of objects is described by a set of rules. In classical P systems, rules enable objects to be created, destroyed or qualitatively changed. In the proposed system, each of these is encapsulated in one type of membranes. In white holes, objects are created, in black holes, they are removed from the system, and at the intersections, they are changed (meaning they can change their state from moving to stopping or vice versa). At the intersections, objects are communicated to other membranes in the system.

Clearly, the behaviour of cars is not strictly deterministic, therefore stochasticity must be introduced at some point. There are two ways to represent the behaviour of a car in the network.

-When an object is introduced to the system at white hole A, a target destination B is chosen from set of available black holes. The shortest path between nodes A and B is calculated and vector of intersection, which must be passed by an object is assigned to it. For each white hole, a set of probabilities is assigned to all available black holes. Therefore, some paths are preferred by the cars then others. -For each intersection, a set of probabilities is assigned to all intersecting roads. When a car is present at the intersection, a road is chosen and car is send to next node of a graph. Each road is associated with a certain probability of being chosen, therefore some links are more preferred by cars than others.

First approach is preferable in case, where the preferences of cars in the network are well known. If we lack the knowledge about drivers preferences and know only vehicle intensities at the roads, second approach seems more suitable. Also, travel time have to be taken into account. Time necessary to travel from one node of a graph to another is represented by a cost of the link. Graph is oriented, therefore travel time between two nodes can differ in both directions.

Last issue is associated with representation of traffic lights. At intersection, only cars approaching from a certain direction are allowed to pass. Within P systems framework, membrane polarizations are suitable tool for selective allowance of objects into a membrane. The polarization of each intersections represents red/green lights and periodically changes, therefore simulating real world situation.

Formal definition of a model

Traffic P system is a construct:

T raf f icΠ = (O, µ, syn, (s (i,j,t) ) (i,j)∈syn;t∈N , R), (2)
where:

-O = {veh, veh s } is an alphabet of objects, veh representing moving vehicle and veh s representing stopping vehicle, µ is a membrane structure, consisting of m membranes, labeled W H, BK, IN T ; representing white holes, black holes and intersections, syn ⊆ {(i, j, t)|i, j ∈ {0, 1, 2, . . . , m}, i = j, t ∈ N} is a subset of synapseslinks between the membranes, where i and j are membranes from the set µ and t is time constant associated with each synapse. -R is a set of rules associated with membranes from µ. The rules differ for membranes labeled as white holes, black holes and intersections.

• White holes -vehicles veh are generated and sent to one of neighboring membranes

. [ ] W H ti -→ [veh W H ] k ,
where {k, W H, t i } ∈ syn and application of a rule takes t i time steps. A label of white hole will be associated with a vehicle veh when approaching intersection k, which will be used to decide, whether a vehicle should stop or pass the intersection. • Black holes -vehicles veh are removed from the system. describes a situation, where stopping car veh s from membrane k is in membrane m with polarization k; it can pass the intersection and continue to membrane n. This is enabled by changing polarizations of membranes.

Global clock for the system are assumed. The polarization of membranes changes according to a given schedule. All other parameters of the system -creation rate at white holes, probabilities and time constants associated with the synapses can also be time-dependent.

XML structure for proposed model

We used XML language for creating topology, creating rules and setting initialization values. The reason why we chose this way is simple, because we can simply and clearly describe all necessary attributes. A code sample is shown in Listing 1.1, where we divide file to 3 parts.

[

  veh] BH t=1 --→ [ ] W H , application of such rule takes 1 time unit.• Intersection -several rules are associated with intersections. Here we describe a situation, where car makes a decision about its target destination at each intersection. Rule:[veh k ] k m ti,pi ---→ [veh m ] ndescribes a situation, where car coming from membrane k approaches an intersection m with polarization k. Since the polarization and label of veh correspond, a veh can pass the intersection and proceed to membrane n, with probability p i . veh will be assigned with label of passed membrane m. Application of rule will take t i time units. Rule:[veh k ] l m t=1,p=1 -----→ [veh s ] l mdescribes a situation, where moving car veh coming from membrane k approaches an intersection m with polarization l. Since the polarization and label of veh do not agree, the vehicle veh must stop and therefore changes to veh s , with probability 1 and taking 1 time step.

Acknowledgement

This work was supported by SGS, V ŠB -Technical University of Ostrava, Czech Republic, under the grant No. SP2012/151 Large graph analysis and processing. This work was also supported by SoftComp project (CZ.1.07/2.3.00/20.0072). SoftComp project is co-nanced by ESF and Czech state budget.

are defined, I/O membranes are chosen etc. 2. Second part describes the initial configuration of membranes. 3. The third part contains the rules for membranes. It consists of left side, right side and additional information about conditions of use, delays, etc. Left side contains items needed for rule to be applied, right side describes the result of rule application.

Listing 1.1. Example XML file <?xml version="1.0" encoding="utf-8"?> <Program Version="1.0"> <Topology> <Membrane ID="1"> <Type Membrane="White"/> <LinkTo ID ="3" ValuePath="5"/> </Membrane> <Membrane ID="4"> <Type Membrane="Black"/> <LinkTo ID ="3" ValuePath="5"/> </Membrane> <Membrane ID="3"> <Type Membrane="Green"/> <PolarizationComing IDs="1,4" Time="10"/> <PolarizationComing IDs="2,5" Time="5"/> <LinkTo ID ="2" ValuePath="10"/> <LinkTo ID ="4" ValuePath="5"/> <LinkTo ID ="5" ValuePath="5"/> </Membrane> </Topology> < Initialization > <Membrane ID="1"> <Item NameID="a" Value="0" /> <Item NameID="b" Value="0" /> </Membrane> </ Initialization > <Rules> <Selected> <Membrane ID="3"> <LeftSide> <Item NameID="a" Value="1" /> </LeftSide> <RightSide> <Item NameID="a" Value="1" /> </RightSide> <Additional> <Probability Value="0.5"/> <DueID Value="4"/> <Polarization Value="true"/> </Additional> </Membrane> </Selected> <Default> <Membrane Typ="White"> <LeftSide> <Item NameID="a" Value="1" /> </LeftSide> <RightSide> <Item NameID="a" Value="1" /> </RightSide> </Membrane> </Default> </Rules> </Program> 7 Concluding Remarks and Future Work A P system for traffic flow simulation was defined and an XML scheme designed to store the configuration of such model was presented. A basic variant of P system was described with additional features -stochasticity, time-dependence, selective membrane permeability -described in more detail.

Presented model was designed to simulate the behaviour of vehicles at the intersections and enables elegant simulation of cars behaviour at the intersections, while ignoring detailed description of behaviour at road segments. This robustness is desirable i.e. for traffic light optimalization.

On the webpages of the authors is possible to find the application, which we develop and use for experiments. (http://wh.cs.vsb.cz/voj189)

In the future, the presented model will be further developed to be able to describe the traffic system in more realistic manner -cars delays, complex intersections and destination preferences are of primal interest.