Ioannis M Stephanakis

Ioannis P Chochliouros
email: ichochliouros@oteresearch.gr

Anastasios Dagiuklas

George C Anastassopoulos

A Novel Rate-Distortion Method in 3D Video Capturing in the Context of High Efficiency Video Coding (HEVC) in Intelligent Communications

Keywords: object oriented coding, H264/AVC, High Efficiency Video Coding HEVC/H.265, higher order motion compensation models, CABAC, 3D television, Rate-Distortion theory

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

Introduction

High Efficiency Video Coding (HEVC) [START_REF] Ohm | Comparison of the Coding Efficiency of Video Coding Standards -Including High Efficiency Video Coding (HEVC)[END_REF][START_REF] Sullivan | Overview of the High Efficiency Video Coding (HEVC) Standard[END_REF] is a video compression standard that is proposed as a successor to H.264/MPEG-4 AVC (Advanced Video Coding) developed by the ISO/IEC Moving Picture Experts Group (MPEG) and the ITU-T Video Coding Experts Group (VCEG) under the name ISO/IEC 23008-2 MPEG-H Part 2 and ITU-T H.HEVC [START_REF]Joint Collaborative Team on Video Coding[END_REF]. MPEG and VCEG have established a Joint Collaborative Team on Video Coding (JCT-VC) in order to develop the HEVC standard. HEVC standard is said to improve video quality, double the data compression ratio compared to H.264, and support 8K UHD and resolutions up to 8192×4320 pixels [4]. Two profiles are supported, namely Main Profile and High Efficiency 10 (HE10). The HEVC Main Profile (MP) is compared in coding efficiency to H.264/MPEG-4 AVC High Profile (HP), MPEG-4 Advanced Simple Profile (ASP), H.263 High Latency Profile (HLP), and H.262/MPEG-2 Main Profile (MP). As of the end of January 2013, ISO/IEC and ITU-T had approved HEVC as ISO/IEC 23008-2 High Efficiency Video Coding and ITU-T Rec. H.265 respectively, as the final draft international standard.

Simulcast as well as Multiview Video Coding are described in the context of the MPEG-4 AVC standard (Version 8: July 2007 (including SVC extension), Version 9: July 2009 (including MVC extension). [START_REF]ITU-T and ISO/IEC JTC 1: Advanced Video Coding for Generic Audiovisual Services[END_REF]. A standard MVC coder consists of N parallelized single view coders. Each of them uses temporal prediction structures and encodes a sequence of successive pictures as intra (I), predictive (P) or bi-predictive (B) frames. Nevertheless MVC is not appropriate for delivering 3-D content for autostereoscopic displays since the bit rate required for coding multiview video with the MVC extension of H.264/AVC increases approximately linearly with the number of coded views. The transmission of 3D video in the Multiview Video plus Depth (MVD) format appears as a promising alternative. ISO/IEC and ITU-T established JCT-3V for the 3D video coding standards as of July 2012. AVC compatible video-plus-depth extension is expected within 2013 [START_REF] Vetro | Overview of the Stereo and Multiview Video Coding Extensions of the H.264/MPEG-4 AVC Standard[END_REF]. The transmission of 3D video in according to Multiview Video plus Depth (MVD) associates encoded depth data representing the basic geometry of the captured video scene with picture frames. Depth data are estimated based on the acquired pictures. They are obtained by the application of depth estimation algorithms and should not be regarded as ground truth. Based on the transmitted video pictures and depth maps, additional views suitable for displaying 3D video content on autostereoscopic displays can be generated using depth image based rendering (DIBR) techniques at the receiver side. All video pictures and depth maps that represent the video scene at the same time instant build an access unit. The access units of the input MVD signal are coded consecutively similar to MVC. The video picture of the so-called independent view inside an access unit is transmitted first directly followed by the associated depth map. Thereafter, the video pictures and depth maps of other views are transmitted. A video picture is always directly followed by the associated depth map. The Coding Tree Unit (CTU) in HEVC replaces the macroblock structure as known from previous video coding standards. It has a size selectable by the encoder and it can be larger than a traditional macroblock. The quadtree syntax of the CTU (see Fig. 1) specifies the size and positions of its luma and chroma Coding Blocks (CB). A Coding Tree Block (CTB) may contain only one Coding Unit (CU) or may be split to form multiple CUs. Thus each CU is characterized by its Largest CU (LCU) size and the hierarchical depth in the LCU that the CU belongs to. It has an associated partitioning into Prediction Units (PUs) and a tree of Transform Units (TUs).

Intra-picture prediction, inter-picture prediction or skip mode are selected at CU level. A PU is basically the elementary unit for prediction and it is defined after the last level of CU splitting. Prediction type and PU splitting type are two concepts that describe the prediction method. Intra prediction allows for symmetric splitting whereas inter prediction allows for both symmetric and asymmetric splitting. At the level of PU, intra-prediction is performed from samples already decoded in adjacent PUs. Such modes as DC (average/flat), angular directions (one from up to 33 as in Fig. 2), Planar Intra Prediction (surface fitting), SDIP (Short Distance Intra Prediction), MDIS (Mode Dependent Intra Smoothing) may be used. Advanced motion prediction (see for example [START_REF] Kordasiewicz | Affine Motion Prediction Based on Translational Motion Vectors[END_REF]) featuring a "merge" mode or the skip mode may be used for inter prediction. Quarter-pixel precision and 7-tap or 8-tap filters are used for interpolation of fractional-sample positions. Multiple reference pictures are used. A deblocking filter similar to the one used in H.264/MPEG-4 AVC is operated in the inter-picture prediction loop. An adaptive loop filter (ALF) is alternatively employed for higher efficiency. The coefficients of ALF are calculated and transmitted on a frame basis and the MMSE estimator is used. For each degraded frame, ALF can be applied to the entire frame or to local areas. Similar transforms as for H.264 (including the discrete sine transform and the Hadamard transform) are used for encoding the residual data. Three different scanning modes (namely zigzag, horizontal and vertical scan) are used to improve the residual coding.

Such novel features as the option to partition a picture into rectangular regions called tiles and wavefront parallel processing (WPP) are introduced into the HEVC standard in order to enhance its parallel processing capabilities.

The core algorithm of Context-Adaptive Binary Arithmetic Coding (CABAC) is based upon that of the H.264/MPEG-4 standard. The number of contexts used in HEVC is substantially less than the number of contexts in H.264/MPEG-4 AVC whereas the entropy coding design allows for better compression.

Motion Compensation Using Higher Order Models and Rate-Distortion Control

Higher order motion models and view synthesis techniques are currently investigated in the literature for sub-pixel motion compensation that may be applied 3D and multiview encoding. Higher order motion models such as the affine model (AMMCP) [START_REF] Wiegand | Affine Multipicture Motion-Compensated Prediction[END_REF], the mesh based MCP [START_REF] Nakaya | Motion Compensation Based on Spatial Transformations[END_REF], the elastic MCP [START_REF] Pickering | Enhanced Motion Compensation Using Elastic Image Registration[END_REF] and View Synthesis by Depth Image Rendering (DIBR) [START_REF] Fehn | A 3D-TV Approach Using Depth-Image-Based Rendering (DIBR)[END_REF][START_REF] Schwarz | 3D Video Coding Using Advanced Prediction, Depth Modeling, and Encoder Control Methods[END_REF] have been proposed as possible extensions of . existing standards. The transformation of one view to another is based on the camera parameters. DIBR uses a view and a depth map for generating arbitrary views. Rough 3D information can be reconstructed based on the depth map so that the transformation due to the disparity effect can be generated. The video compression with DIBR becomes the compression of depth image instead of the correlations between views. The elastic MCP model estimates the spatial transformation parameters between the predicted block with translational MV and the current block.

MPEG-4 HEVC

The motion parameters are encoded in the bit stream. The pixel location of the predicted block (x i , y i ) is transformed to (x' i , y' i ) by the following equations,

    2 / 1
) , ( '

P l i i l l i i y x m x x  and      P P l i i l l i i y x m y y 1 2 / ) , ( '  ( 1
)
where P is the number of elastic parameters used, m l is the elastic parameter and ) , (

i i l y x 
is the basis function. The elastic MCP model uses discrete cosines as the basis functions and it is defined by

                   y blocksize v y x blocksize u x y x y x i i i i P l i i l 2 ) 1 2 ( cos 2 ) 1 2 ( cos ) , ( ) , ( 2 /     , ( 2
)
where 1 2 /    v u P l and 1 2 / ) , ( 0    P v u
. Conventional motion vector is equal to motion parameters [m 1 m P/2+1 ] in such a case.

The minimization of a Lagrangian cost function for motion estimation was proposed in [START_REF] Sullivan | Rate-Distortion Optimized Motion Compensation for Video Compression Using Fixed or Variable Size Blocks[END_REF] (see Fig. 3). Given a reference picture list R and a candidate set M of motion vectors, the motion parameters for a block s k , which consists of a displacement or motion vector m=[m x , m y ] and, if applicable, a reference index r, determine the coding mode for coding a block of samples, such as a macroblock or a coding unit. Additional features of a coding mode are the intra or inter prediction modes or partitions for motion-compensated prediction or transform coding including the quantization step. Given the set of applicable coding modes for a block of samples s k , the used coding mode is chosen according to where the distortion term D k (c) represents the SSD between the original block s k and its reconstruction s' k , that is obtained by coding the block s k with the mode c. The term R k (c) represents the number of bits (or an estimated thereof) that are required for representing the block s k using the coding c for the given bitstream syntax. It includes the bits required for signaling the coding mode and the associated side information (e.g. motion vectors, reference indices, intra prediction modes and coding modes for sub-blocks of s k ) as well as the bits required for transmitting the transform coefficient levels representing the residual signal. A coding mode is often associated with additional parameters such as coding modes for sub-blocks, motion parameters and transform coefficient levels. While coding modes for sub-blocks are determined in advance, motion parameters and transform coefficient levels are chosen according to Eq. 3. For calculating the distortion and rate terms for the different coding modes, decisions for already coded blocks of samples are taken into account (e.g. by considering the correct predictors or context models). 4. The design of binarization schemes relies on a few basic code trees, whose structure enables a simple on-line computation of all code-words without the need for storing tables. There are four such basic types [START_REF] Marpe | Context-Based Adaptive Binary Arithmetic Coding in the H.264 / AVC Video Compression Standard[END_REF] namely the unary code, the truncated unary code, the k-th order Exp-Golomb code and the fixed-length code. There are binirization schemes based on a concatenation of these elementary types. The nodes located in the vicinity of the root node of a binary tree are the natural candidates for being modelled individually, whrereas a joint model should be assigned to all nodes on deeper tree levels corrresponding to the tail of the probability density function. According to CABAC, given a predefined set T of past symbols, a so-called context template, and a related set C={0,…,C-1} of contexts, contexts are specified by a modeling function C T F  : operating on the template T. For each symbol x to be coded, a conditional probability p(x|F(z)) is estimated by switching between different probability models according to the already coded neighboring symbols T z  . After encoding x using the estimated conditional probability p(x|F(z)), the probability model is updated with the value of the encoded symbol x. Thus p(x|F(z)) is estimated on the fly by tracking the actual source statistics. The entity of probability models used in CABAC can be arranged in a linear fashion such that each model can be identified by a unique so-called context index .

)) ( ) ( ( min arg * c R c D c k k C c k      (3)

A Scalable Context-based Adaptive Model for Encoding Motion Compensation Parameters in Multi-view 3D Systems

A Context-based Adaptive Model for Correlated Motion Compensation Parameters in Multiple Views

The proposed encoding scheme investigates context-based adaptive models for correlated multilinear object entities between multiple views [START_REF] Stephanakis | A Multiplicative Multi-linear Model for Inter-Camera Prediction In Free View 3D Systems[END_REF]. It is a model that does not depend upon the parameters of the recording cameras. It is applicable to such object entities as transform elements and/or motion vector parameters and it scales in a straightforward fashion from stereo to multiple free views. It combines adaptivity regarding the context templates with updates of the correlated orthonormal features of the object entities through all views per GOP. Let us assume an N th -order tensor A, which resides in the tensor multi-linear space

N I I I R R R     2 1
where

N I I I R R R , , 2 1



are the N vector linear spaces of multiview system featuring N views. The "k-mode vectors" of A are defined as the I k -dimensional vectors obtained from A by varying its index in k-mode while keeping all other indices fixed [START_REF] De Lathauwer | A Multilinear Singular Value Decomposition[END_REF][START_REF] Lu | MPCA: Multilinear Principal Component Analysis of Tensor Objects[END_REF]. Multiview motion prediction using cross-view prediction vectors may be defined per GOP through a video object that we call Motion Prediction video-Object and is denoted as MPO. Let us define a group of motion vector parameters pertaining to k-view within some GOP as

  k k I k k k M M M ) ( 2
           ) ( 1 1 2 / 2 / 1 ) , , ( ] [ ] [ ) , , ( ) , , ( k I
i k i block block k i y x P P l l k l P l l k l block block block block k M t n m v v m m t n m t n m    v v . ( 4
)
Motion Prediction video-Object (MPO) holds the orthonormal cross-view prediction vectors. It is defined as,

) ( ) 3 ( 3 ) 2 ( 2 ) 1 ( 1 2 1 2 1 1 1 1 ) ( ) ,...; , ( ) , , ( 2 1 1 1 2 2 N N N i i i I i I i I i block block t S M M M t i i s t n m MPO N N N M M M M                 (5)
Unfolding MPO along the k-mode is denoted as

) ( ) ( 1 1 1 Ν k k k I I I I I k R x x x x x x      MPO , (6)
where the column vectors of MPO (k) are the k-mode vectors of MPO.

Unfolding the Multi-view Video Plane of an N-view system along the k-mode view results the following matrix representation

  T ) 1 ( ) 1 ( ) ( ) 2 ( ) 1 ( ) ( ) ( ) ( ) ( ) (           k N k k k k k t S t M M M M M M MPO   . . , (7)
where  denotes the Kronecker product. The core tensor S (in a representation similar to the one described in Eq. 8) is analogous to the diagonal singular value matrix of the traditional SVD and coordinates the interaction of matrices to produce the original tensor. Matrices M (k) are orthonormal and their columns span the space of the corresponding flattened tensor denoted as MPO (k) . The objective of MPC analysis for predetermined dimensionality reduction is the estimation the N projection matrices

{ N k R t k k P I k ,..., 1 , ) ( ~) (   x M
} that maximize the total tensor scatter [START_REF] Lu | A Survey of Multilinear Subspace Learning for Tensor Data[END_REF],

  2 ... , ) ( ) ( ) ( max arg ,..., 1 ), ( ~) ( (2) (1) t Δ t MPO N k t N k    M M M M , (8)
where

c P N k k   1
. One may estimate adaptively the cross products between motion vector parameters of different video views by maximizing sequentially Eq. 8 for all modes. Its decomposition into non-negative parts and the application of a multiplicative update rule that maintains orthonormality [START_REF] Yang | Multiplicative Updates for Non-negative Projections[END_REF][START_REF] Zhang | Effective Multiplicative Updates for Non-negative Discriminative Learning in Multimodal Dimensionality Reduction[END_REF] is suggested in [START_REF] Stephanakis | A Multiplicative Multi-linear Model for Inter-Camera Prediction In Free View 3D Systems[END_REF].

The following algorithm outlines in detail the steps of the proposed encoding approach for multi-view systems:

Select reference frames for GOP and initialize motion vector parameters per GOP I -Estimate the average motion vectors [v x v y ] that are common for all views according to Eq. 4. II -Subtract average values and estimate the sets of motion vector parameters that are strongly correlated. Carry out prediction using adjacent macroblocks and solve Eq. 8 iteratively (on the fly) for optimal MPOs using multiplicative updates (or operational updates). Assume (Fig. 5) that

) ( ) ( ' ) , ; ( ) ( ) ( ) ( ) , ( ) ( ~) ( ) ( ) ( )
w b M w b M Φ M            ( 9
)
where t is the frame index, w(t) stands for a white process featuring zero mean and (t)denotes some type of operation at t like element-byelement multiplication, component rotation, permutation etc.

III-

Numerical Simulations

Numerical simulations for the proposed encoding method have been carried out for the image sequences used for video view interpolation as described in [START_REF] Zitnick | High-quality video view interpolation using a layered representation[END_REF]. Each sequence is 100 frames long. The camera resolution is 1024x768 and the capture rate is 15fps. The frames of the ballet sequence in [START_REF] Zitnick | High-quality video view interpolation using a layered representation[END_REF] are used to obtain the numerical values presented in this section. The multi-view GOP for the numerical simulations consists of initial four (4) frames. A fixed size macroblock featuring dimensions of 8x8 pixels is used to carry out motion compensation. The average of the first frame of View 1 and of the first frame of View 2 is used as a reference. The elastic MCP model as described in Eq. 2 uses sixteen (16) discrete cosines as the basis functions (parameter P equals 32). The context and the distributions for the first stereo GOP of the sequence are illustrated in Fig. 6 Entropy encoding is used based upon the probabilities depicted in Fig. 6.a. The 32x32 sigma values ordered according to magnitude as obtained from the proposed algorithm are given in Fig. 6.b. The histogram of the highest sigma values over all macroblocks is presented in Fig. 6.c. The magnitude of the residual frames corresponding to View 1 and View 2 are marginally decreased as linear correlated components between views are taken into account. Peak-Signal-to-Noise-Ratio (PSNR) values corresponding to the two views are 31.9 dB and 31.0 dB respectively when the two (c=2) most correlated linear components are used. The additional overhead is estimated to 0.11 bits/pixel. Original and residual frames are given in Fig. 7. PSNR values are slightly improved when the six (c=6) most correlated linear components are taken into account. PSNR value for View 1 equals 31.9 dB whereas PSNR value for View 2 equals 31.2 dB for an estimated additional overhead of 0.17 bits/pixel. Nevertheless savings on the bits/pixel required for underlying motion compensation could be made should the proposed approach be applied.

Conclusion

A novel scalable approach to multi-view video encoding based on the so-called MPO structure is proposed. It takes advantage of the correlated linear components between views and it requires no prior knowledge of the capturing cameras in space. Strongly correlated linear components, i.e.

Fig. 1 .Fig. 2 .

 12 Fig. 1. Structure of Coding Tree Units (CTU). Prediction and Transform Units (PU and TU) are depicted Fig. 2. Intra-picture prediction for boundary samples of adjacent blocks (MPEG-4 and HEVC)

Fig. 3 .Fig. 4 .Fig. 5 .

 345 Fig. 3. Rate control model

Fig. 6 Fig. 7

 67 Fig. 7 Views 1 and 2 and residual frames for joint motion vector compensation (the two largest s(i 1 , i 2 ,…;t) values have been used) multiplicative per GOP (or frame or slice). The correlated linear components are assigned initial values per GOP (or frame or slice) whereas subsequent operations (updates and permutations) are defined upon index position. Context templates corresponding to different distributions of the sigma correlation parameters may be selected according to the CABAC model. Updating adaptively indexed orthonormal basis functions in conjunction with CABAC increases encoding efficiency and allows for the incorporation of SVD and MPC transforms into the existing approaches. Initial numerical results indicate that the proposed method yields improvements in encoding multiple views in MPEG standards under development like the HEVC.

  It improves the compression capabilities for dependent video views and depth data. As of the end of January 2013, ISO/IEC and ITU-T had approved HEVC as ISO/IEC 23008-2 High Efficiency Video Coding and ITU-T Rec. H.265 respectively, as final draft international standard.

	 Encoder Control
	 View synthesis algorithms
	Final draft amendment (FDAM) for HEVC for stereo and 3D is expected in 2015.
	Advanced Motion Vector
	Prediction (AMVP) is used. Work on HEVC 3D and scalable extensions has been
	currently under development. It focuses on
	 Coding of Independent Views (2D video coding)
	 Coding of Dependent Views (DCP, inter-view motion/residual prediction)
	 Coding of Depth Maps

2 Higher Order Models for Motion Compensation and Encoder Control 2.1 Novel Video Encoding Features of the HEVC/H. 265 Standard

	at each "clean" random access point
	into a video sequence) is coded using only intra-picture prediction. Transform
	coefficients from region-to-region within the same picture are spatially predicted but
	there are no dependencies upon other pictures. The remaining pictures of a sequence
	or frames between random access points are encoded using temporally-predictive
	coding modes. The encoding process for inter-picture prediction consists of choosing
	motion data comprising the selected reference picture and motion vectors (MV) that
	are applied for predicting the samples of each block. Residual frames are transformed
	by a linear spatial transform. The transform coefficients are scaled, quantized, entropy
	coded and transmitted together with the prediction information. Novel features of the
	HEVC standard are outlined as follows:
	 Larger and more flexible coding, prediction, and transform units.
	 Improved mechanisms to support parallel encoding & decoding.
	 More flexible temporal prediction and scanning structure.
	 More accurate intra prediction approach and directions/modes.
	 More accurate motion parameters (including merge mode) and sub-pixel
	prediction.
	 Inclusion of non-square transform and allowing asymmetric motion
	prediction.
	 More flexible transform, choice of DST, and no-transform option.
	 Rate-distortion optimized quantization (RDOQ).
	 Improved in-loop filters, including the new sample adaptive offset

The video coding layer of HEVC employs the same "hybrid" approach used in all video compression standards since H.261. Each picture is split into block-shaped regions, with the exact block partitioning being conveyed to the decoder. The first picture of a video sequence (and the first picture

  -is carried out once with respect to one of the views or a linear combination of a selected subset. The motion vectors pertaining to k-view according to the elastic model may be decomposed as follows,

		M	(	)			1									. Prediction of motion vectors -which
	are denoted as	v	(	block m	,	block n	,	t	)		[ v	x	v	y	]	for a block indexed by	(	, block n block m	)	at
	t

  (

	k	t	t	View	k	t	t	t	k	t	t	t	F View t	k	k	t	t	t	t	t

Acknowledgments. This research work has been funded by the LiveCity

European Research Project supported by the Commission of the EC -Information Society and Media Directorate General (FP7-ICT-PSP, GA No.297291).


Salvador Valencia Rodríguez

Herna Lydia Viktor
email: hviktor@uottawa.ca

A Personalized Location Aware Multi-Criteria Recommender System Based on Context-Aware User Preference Models

Keywords: Recommender Systems, Location Aware, Multi-Criteria, Preference Models, Personalization

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

Introduction

Consumers often find themselves in situations where they have to choose one item over others. For instance, they may wish to decide which movie to view, which book to read, what items to buy, and so forth. However, due to the advance in technology and the large amount of information available in databases, our options have dramatically increased. We have now reached a point where we have thousands of options at our fingertips, through the use of web-based systems. Consequently, in order to address this information overload, and aid consumers through these daily decisionprocesses, Recommender Systems have been developed. These systems aim to provide personalized services to each user, showing them only the information that they are most likely to be interested into [START_REF] Ricci | Recommender Systems Handbook[END_REF].

A number of techniques to predict the best items have been proposed. While some of these Recommender System implementations have been successful in many domains, a number of challenges still remain. Most implementations consider only single criteria ratings, and consequently are unable to identify why a user prefers an item over others. Some systems classify the user into one single group or cluster, an approach that has limited use, since real world users share commonalities in different degrees with diverse types of users. Finally, other systems require a large amount of previously gathered data about users' interactions and preferences, in order to be successfully applied.

This work introduces an algorithm to overcome these previously mentioned disadvantages. Our algorithm, as presented in this paper, builds user preference models considering multiple criteria. That is, we include the users' special needs, and context in the decision making process. This enables our system to clearly identify most important criteria, for a specific user when selecting an item over the others, and correspondingly create accurate recommendations that match his/her preference's value scale. Moreover, by including the user context as part of the recommendation process, the system is able to produce different types of recommendations to the same user, depending on his/her current context. Additionally, our algorithm exploits the information contained in a large-scale demographic database to generalize the information as provided by the user. This is done by clustering the user into one or more demographic groups. This aspect allows our system to leverage commonalities between similar user types, and create richer user profiles without the need of previously stored data about other users' interactions.

This paper is organized as follows. Section 2 details our recommender system algorithm. Section 3 describes our case study and Section 4 concludes.

Personalized Location-Aware System

We introduce a Personalized Multi-Criteria Context-Aware Recommender-System, which has been designed to achieve the following goals. Firstly, we aim to identify the user preferences, as based on multiple criteria that lead him/her to select an item over others, and correspondingly potentially generate more accurate recommendations. Further, our objective is to model and consider the user context during the recommendation process, instead of using it only as a filter. That is, we follow a contextaware approach, where the context is used as one of multiple criteria for the creation of preference models [START_REF] Lee | Location-Aware Agent Using Data Mining for the Distributed Location-Based Services[END_REF][START_REF] Schilit | Disseminating active map information to mobile hosts[END_REF]. Finally, our ultimate goal is to produce accurate recommendations by generalizing user profiles without gathering data from previous users' interactions. We use a demographic generalization technique to exploit the information contained in large-scale demographic databases that encompass interests and preferences of similar people [START_REF] Krulwich | Lifestyle Finder. Intelligent User Profiling Using Large-Scale Demographic Data[END_REF].

Our algorithm is based on the creation of a multiple criteria user preference model considering the following four criteria. These are 1) the probability that a user prefers an item, given the probabilities that the user belongs to each of the demographic clusters identified in a large-scale demographic database; 2) the weight of each item, based on its attributes and the weights the user assigns for each of them; 3) the distance between the user's current location and the item's location; and, 4) the time the user would require to reach each item, based on his/her mobility level, together with other geo-spatial constraints (such as a river) and average route times. The user preference model is then used to assess the utility of each item for the user. Finally, the system will recommend to the user the items with the highest estimated utility. The following figure presents the high-level schematic overview of our algorithm.

Create the User Profile

We create a user profile consisting of two dimensions. The first is the set of all possible demographic clusters where the user potentially belongs to. The second aspect concerns how much importance a specific user assigns to an item. These dimensions may be expressed as the following two vectors. DemogClusters(u), which represents the possible demographic clusters where a user may belong, along with the probability of belonging to each of them. The second is AttributeWeights(u), which contains the possible items' attributes along with the user-assigned weight for each of them. The set of demographic clusters where a user may belong (k clusters), that form the DemogClusters(u) vector is obtained by applying a demographic generalization technique. First, the algorithm select the m demographic variables that best differentiate each cluster from the others using feature selection techniques. Second, the user is queried to obtain his/her information for the selected demographic variables. Third, our system obtains the probability of a user belonging to each demographic cluster, based on the user given information and the users' distribution (included in the database) for each possible value for each demographic variable. Finally, the set of possible demographic clusters where the user may belong, are those in which the product of the users' distributions for the user given information is higher than a selected threshold. The following equation presents how these clusters are obtained.

  , 1 : ( ( , )) 1 x m x u C j j j PossibleClusters u C P DV Val DV u Threshold x N                  (1)
Here, DV stands for Demographic Variable;

, represents the value of the demographic variable j given by user u; and finally, , % expresses that there is a Z% probability (obtain from the users' distributions) that user u belongs to cluster x, given that the value of the demographic variable j for that user is Y.

Once the demographic clusters a user could belong to have been selected, the probabilities already obtained from equation ( 1) may be considered as the probabilities of belonging to each of them. However, this would imply classifying the user based only on his/her demographic information, and not considering other data that might also be included the database (e.g. information about leisure, shopping, media preferences, hobbies, etc). Therefore, for databases where this type of information is available, we propose the following technique to obtain the probabilities of belonging to each of the selected demographic clusters. First, we obtain the Cartesian product of the nondemographic information included in the database for each of the previously selected clusters (e.g. Leisure Х Shopping X Media). Second, we apply the k-means data mining clustering algorithm over the Cartesian product dataset, to learn the demographic clusters and create the same number of groups as previously selected clusters (k groups).

The results of the k-means algorithm are k groups, each labeled by the most representative non-demographic information it contains. The k-means algorithm also creates a table, with the distribution of instances from each demographic cluster classified into the newly created groups, as shown in Table 1.

Table 1. Instances from demographic clusters classified in the newly created groups

In Table 1, Inst(Cx) represents the total number of instances belonging to the demographic cluster x, Inst(Gy) represents the total number of instances assigned to the newly created group y, and Inst(Cx,Gy) represents the number of instances from the demographic cluster x that were assigned to the newly created group y.

Third, we present to the user the obtained newly created groups, so he/she can select the one (or more) group(s) that best define him/her (p selected groups), based on the non-demographic information that labels each of them. Fourth, we obtain probabilities for each of the k demographic clusters, using the following equation:

      1 1 If , 0 : 0 , 1 , *100 1 Otherwise : ( ) p x j j x p x j j j Inst C G P u C x k Inst C G p I n s tG                         (2)
Once we have obtained the possible demographic clusters where a user may belong and computed the probability for each of them, the DemogClusters(u) vector is complete. Equation [START_REF] Siskos | UTA Methods[END_REF] shows the structure of this vector, formed by pairs of the form , , , where represents the demographic cluster x, and , the probability that the user u belongs to cluster x.

          1 1 2 2 , , , , , , , , , k k DemogClusters u
C P u C C P u C C P u C         (3)
The second vector included in the user profile (AttributeWeights(u)) consists of pairs of the form , , , where represents the item's attribute x, and , is the weight that user u assigns to attribute x. Equation ( 4) shows the structure of this vector, where r represents the selected number of items' attributes that best differentiate the items between them:

          1 1 2 2 , , , , , , , , , r r

AttributeWeights u A W u A A W u A

A W u A         (4)

Location Awareness and Mobility Level

The second phase of the algorithm consists of obtaining the user's mobility level and current location. Consequently, the user is asked to choose the means of transport (i.e. mobility level) that best describes his/her current situation (e.g. driving a car, riding a bike, wheelchair, or walking). Based on the user selection, the routes to be considered when predicting the distance and time to each item are selected. Subsequently the user's current latitude and longitude spatial coordinates are obtained from Google Maps®, through its application programming interface (API), using the current user address.

Assessing the Items' Utilities Considering the User Profile

The third phase of our algorithm is divided into two stages. The first stage entails obtaining an utility for each item (Demographic Utility (DU)), based on the previously created vector DemogClusters(u). However, in order to do so, all the items in the database need to be mapped to one or more demographic clusters. Consequently, we assign each item to one or more categories, which in turn are linked to one or more demographic clusters. Once all the items are related to one or more demographic clusters through one or more categories, we compute the probability that a user might prefer each of these categories. As shown in equation ( 5), these probabilities are calculated based on the probabilities of belonging to each of the k previously selected demographic cluster.

 

1 If : ( , ) , Otherwise : 0

k x i i x i ItemCat C P u C P u ItemCat        (5)
Here, , represents the probability that user u prefers the item category x. Subsequently, these probabilities are used to obtain the demographic utility (DU) for each item. This DU is computed as the average of the probabilities that the user prefers each category where the item belongs to. The DU of an item is obtained by using equation [START_REF] Lee | Location-Aware Agent Using Data Mining for the Distributed Location-Based Services[END_REF], where q represents the number of identified items' categories.

  1 If : ( , ) 1 ,
Otherwise :0

q x i i x i x I ItemCat P u ItemCat DU I u Numof Cat of I        (6)
The second stage involves assessing a utility for each item (Weighted Utility (WU)), considering now the previously created AttributeWeights(u) vector. The WU for each item is calculated as the sum of each user given weight for all the attributes of an item, as shown in equation [START_REF] Schilit | Disseminating active map information to mobile hosts[END_REF].

  1 If : ( , ) , Otherwise : 0

r i x i x i A I W u A WU I u        (7)

Assessing the Items' Utilities Considering the User Location

The fourth step in the algorithm involves assessing a utility for each item, based on the user's current location and mobility level. In this phase, each item is assessed under two additional criteria: the distance from the user (Distance Utility (DisU)), and the time it would take the user to reach the item (Time Utility (TU)), considering the user's mobility level. It is important to note that these new utilities, expressed in kilometers and minutes from the user, are inversely proportional to the item utility for the user. That is, the closer an item is the higher its utility for the user. Finally it is noteworthy to mention that while the TU is obtained from Google Maps®; the DisU is obtained by means of SQL spatial queries, performed over a spatial database using the user coordinates.

Creating the User Preference Model

The utilities that were obtained during the third and fourth phases represent the four criteria considered to create the user preference model, which in turn represents the importance to the user. To obtain the user preference model we apply the UTA* algo-rithm. This algorithm is a regression-based technique that infers preference models from given global preferences (e.g. previous user choices) [START_REF] Lakiotaki | Multi-Criteria User Modeling in Recommender Systems[END_REF][START_REF] Siskos | UTA Methods[END_REF]. For more information on this algorithm, we recommend [START_REF] Lakiotaki | Multi-Criteria User Modeling in Recommender Systems[END_REF][START_REF] Siskos | UTA Methods[END_REF][START_REF] Siskos | UTASTAR: An ordinal regression method for building additive value functions[END_REF][START_REF] Jacquet-Lagrèze | Assessing a set of additive utility functions for multicriteria decision making: The UTA method[END_REF].

In order to apply this algorithm we first need to obtain the user's weak preference order, which represents the user preferences for some items, after considering the four criterions for each of them. Consequently, the system randomly selects ten items, along with their four utilities, to present to the user, who is then asked to rate them in a non-descending order. The result of the UTA* algorithm is a vector of four values, which represent the weights of each of these four criterions for that user. These weights are calculated using linear programming techniques. The resulting vector from the UTA* algorithm has the following form:

  { , , , } DU WU DisU TU PM u W W W W  (8)

Recommending the Items that Best Match the User Preferences

The last phase of the algorithm involves assessing an integrated, final utility for each item, using the previously obtained user preference model. It follows that, to apply the previously obtained weights from the model, the values of the four different utilities must be normalized so they are in the same range {0,…,1}, otherwise the final utility would be biased toward the criterion with the highest scale. Therefore, each utility is divided by its highest value among all the items. Moreover, since the distance and time utilities are inversely proportional to the user utility, after dividing their values by their highest utility the resulting value is subtracted from 1, and the remainder is the utility to be considered. The final integrated item's utilities are obtained by applying the following equation:

            , , , ( , ) , 1 1 max( ( )) ( ) max( ) max( ( )) x x x x x D U W U D i s U T U DU I u WU I u DisU I u TU I u U I u W W W W DU u max WU u DisU u TU u                     (9)

Experimental Evaluation

As a case study, we created PeRS, a Personal Recommender System that recommends events to attend to consumers. The source datasets used in this case study are the following. The Ottawa Open Data Events database is our items dataset, which contains the events that will take place within the National Capital Region of Canada. By events we refer to festivals/fairs, film/new media, galleries, music, theater, and so on [START_REF]Ottawa Open Data[END_REF]. Our large-scale demographic database is the PRIZM C2 Database that classifies Canada's neighborhoods into 66 unique lifestyle types, providing insights into the behavior and mindsets of residents [10]. The Ottawa Open Data Geospatial databases contain spatial information to geographically locate different elements of interest (i.e. rivers, buildings, museums, municipalities, wards, roads, cities and country areas) [START_REF] Berenson | Basic Business Statistics: Concepts and Applications[END_REF].

Experimental Design

We evaluated our algorithm using an offline experiment. This type of experiment was chosen because it focuses on the recommendation technique that is being used, rather than the system interface, which makes it highly suitable for our environment [START_REF] Ricci | Recommender Systems Handbook[END_REF].

Since the design and implementation of a Recommender System depends on the specific requirements, goals, and sources of information, the properties used to evaluate a system should be selected according to each application domain. The selected properties to be evaluated in our offline experiment were the system's prediction accuracy, item and user coverage, and confidence. (For more information on other properties to evaluate and compare Recommender Systems, we recommend [START_REF] Ricci | Recommender Systems Handbook[END_REF].)

We performed a pilot study with 30 human subjects to obtain, according to the central limit theorem, a punctual estimation of the standard deviation of the population [START_REF] Berenson | Basic Business Statistics: Concepts and Applications[END_REF]. Table 2 presents the structure of the questionnaire answered by the human subjects that took part in the experiment. This questionnaire was designed to gather the required data to create the users' preference models, and to evaluate the recommendations.

Table 2. Structure of the questionnaire applied to the users

Section

Gathered Information Personal Information

The user is asked 9 personal questions, along with his/her mobility level Preferences

The user is asked to provide weights to the identified events' attributes Activities

The user is asked to select the groups of activities where he/she identifies the most.

Events

The user is asked to rate 10 randomly selected events, considering the information for all the criteria.

Experimental Results

Table 3 presents a summary of the results obtained. Our analysis is centered on four aspects. Firstly, our goal was to analyze the distribution of the subjects considered in the experiment. Secondly, we aimed to evaluate the system accuracy and the ability to cluster the subjects into the most adequate demographic groups. Thirdly, we learned the system accuracy from the user profiles. Finally, we evaluated the system coverage, in terms of user and item space coverage.

From the results, as shown in Table 3, we may conclude that, based on their demographic information, 93% of the subjects that partook in experiment have differences in their profiles. This allowed us to test the system for a wide-spread distribution of users, and therefore potentially draw more reliable conclusions for the general population. Additionally, we concluded that the system was able to narrow down the possible matching clusters where a subject may belong, using only the selected demographic variables asked to the subjects. Furthermore, the subjects identified themselves within 86% of the clusters selected by the system, considering the characteristic activities of the people within each group. This indicates that the system accurately classified the subjects into the correct demographic groups.

The system has a predictive accuracy between 75% and 82% with a confidence of 95%. This accuracy represents the effectiveness of the system to identify the user preferences and rank the items in the same way he/she would approach this task. Moreover, for 80% of the subjects considered in our experiment the system was able to rank the items with at least 70% match with the way they would have proceeded (i.e. user space coverage). In addition, the system recommended 81% of the items in the database to at least one user (i.e. item space coverage).

Finally, from the results obtained from the performed classification and statistical analyzes, we concluded that the accuracy of the produced recommendations does not solely depend on the user profile. Therefore, it could potentially be used for a wide range of the population and produce equally accurate recommendations. We will explore the validity of this observation in our future work. The users present an equal distribution in most of the considered personal information. Correspondingly this wide spread distribution of users allows the system to be tested for different user types, and therefore draw reliable conclusions for a wide range of the population.

Table 3. Experimental results

Property

Measured Value Interpretation

Narrowing Percentage

The system identified an average of 4 out of 66 possible clusters where a user may belong

The system is able to narrow in a 95.38% the possible demographic groups where a user may belong, during the demographic generalization process. 86.66% of the users identified themselves with the activities included in the demographic groups selected by the system.

Accuracy

Prediction Accuracy and

{75.36 % -82.96%}

Confidence: 95%

The system recommends items that match the user preferences with an average of 79.16%, for 95% of the population.

Learning the

System Accuracy

Based on the results obtained from applying classification and statistical regression analyzes, using a decision tree, we conclude that the system accuracy doesn't depend on the user profile. Therefore it follows that the system can produce equally accurate recommendations for a wide range of the population.

User Space

Coverage

For 80% of the users the system created models that can rank the events with at least a 70% match with how the users would have ranked them.

Item Space

Coverage 81.97% of the items in the database were recommended to a user at least once, with at least 50% match of his/her preferences.

Conclusion

This paper presents a Recommender System that utilizes an algorithm to generate user models which includes the user context and preferences in the decision model. Our

Location Aware Multi-Criteria Recommender Systems 10 algorithm, as shown in the experiments, accurately identifies the demographic groups where a user may belong. Our results indicate that we are able to produce highly accurate recommendations for a wide range of the population. Further, we are able to construct accurate user profiles. Our algorithm has the following noteworthy characteristics. It uses more than one decision criteria and is therefore able to identify what it is most important for a user when selecting an item over the others. It creates user preference models, consequently the system is capable to understand the user interests and preferences, and is not biased from previous user actions. It generalizes the user profile, thus the system is capable to reason from the commonalities between demographic groups, and it doesn't require previous gathered information about other users' interactions. Importantly, it is location aware and as a result the system varies the recommendations as the user location or level of mobility changes. Our future work will focus on producing accurate preference models in an imbalanced preference setting. This would imply that our algorithm should distinguish user preferences, even when they may be very similar. Additionally, it would be interesting to evaluate the prediction accuracy of our system considering a larger number of users and different locations. Finally, it would be worthwhile to compare the obtained results from our system, against other linear and non-linear recommendation techniques [START_REF] Klašnja-Milićević | The use of Nonlinear Manifold Learning in Recommender Systems[END_REF].

Fig. 1 .

 1 Fig. 1. High-level schematic overview of our algorithm

Distribution 93 .

 93 33% of the users in the experiment have different profiles (based on their personal information).

  30 users identified themselves, in every one of the selected demographic groups for him/her.


Georgios Sermpinis
email: georgios.sermpinis@glasgow.ac.uk

Charalampos Stasinakis
email: c.stasinakis.1@research.gla.ac.uk

Andreas Karathanasopoulos
email: a.karatahnasopoulos@uel.ac.uk

Kalman Filter and SVR combinations in forecasting US Unemployment

Keywords: Forecast Combinations, Kalman Filter, Support Vector Regression, Unemployment

The motivation for this paper is to investigate the efficiency of a Neural Network (NN) architecture, the Psi Sigma Network (PSN), in forecasting US unemployment and compare the utility of Kalman Filter and Support Vector Regression (SVR) in combining NN forecasts. An Autoregressive Moving Average model (ARMA) and two different NN architectures, a Multi-Layer Perceptron (MLP) and a Recurrent Network (RNN), are used as benchmarks. The statistical performance of our models is estimated throughout the period of 1972-2012, using the last seven years for out-of-sample testing. The results show that the PSN statistically outperforms all models' individual performances. Both forecast combination approaches improve the statistical accuracy, but SVR outperforms substantially the Kalman Filter.

Introduction

Many applications in the macroeconomic literature aim to derive and compare information from econometric models' forecasts. For that reason, forecasting competitions of linear and non-linear architectures are common and focus on numerous time series, such as unemployment, inflation, industrial production, gross domestic product etc. The artificial NNs are computation models that researchers include in such macroeconomic forecasting schemes, because they embody promising data-adaptive learning and clustering abilities. The motivation for this paper is to investigate the efficiency of a Neural Network (NN) architecture, the Psi Sigma Network (PSN), in forecasting US unemployment and compare the utility of Kalman Filter and Support Vector Regression (SVR) in combining NN forecasts. An Autoregressive Moving Average model (ARMA) and two different NN architectures, a Multi-Layer Perceptron (MLP) and a Recurrent Network (RNN), are used as benchmarks. The statistical performance of our models is estimated throughout the period of 1972-2012, using the last seven years for out-ofsample testing. The results show that the PSN statistically outperforms all models' individual performances. Both forecast combination approaches improve the statistical accuracy, but SVR is substantially better than the Kalman Filter.

Section 2 is a short literature review and Section 3 follows with the description of the dataset used in this application. Sections 4 and 5 give an overview of the forecasting models and the forecast combination methods implemented respectively. The statistical performance of our models is presented in Section 6. Finally, some concluding remarks are summarized in Section 7.

2

Literature Review

Forecasting unemployment rates is a very popular and well documented case study in the literature (see amongst others Rothman [START_REF] Rothman | Forecasting Asymmetric Unemployment Rates[END_REF], Montgomery et al. [START_REF] Montgomery | Forecasting the U.S. Unemployment Rate[END_REF] and Koop and Potter [START_REF] Koop | Dynamic Asymmetries in U.S. Unemployment[END_REF]). Swanson and White [START_REF] Swanson | A Model Selection Approach to Real-Time Macroeconomic Forecasting Using Linear Models and Artificial Neural Networks[END_REF] forecast several macroeconomic time series, including US unemployment, with linear models and NNs. In their approach, NN architectures present promising empirical evidence against the linear VAR models. Johnes [START_REF] Johnes | Forecasting unemployment[END_REF] reports the results of another forecasting competition between linear autoregressive, GARCH, threshold autoregressive and NNs models, applied to the UK monthly unemployment rate series. In his application, NNs are superior when the forecasting horizon is 18 and 24 months ahead, but fail to outperform the other models in shorter forecasting horizons.

Liang [START_REF] Liang | Bayesian neural networks for nonlinear time series forecasting[END_REF] applies Bayesian NNs in forecasting unemployment in West Germany and shows that they present significantly better forecasts than traditional autoregressive models. Teräsvirta et al. [START_REF] Teräsvirta | Linear models, smooth transition autoregressions, and neural networks for forecasting macroeconomic time series: A reexamination[END_REF] examine the forecast accuracy of linear autoregressive, smooth transition autoregressive and NN models for 47 monthly macroeconomic variables, including unemployment rates, of the G7 economies. The empirical results of their study point out the risk for implausible NN forecasts at long forecasting horizons. Nonetheless, their forecasting ability is much improved when they are combined with autoregressive models. This idea of combining forecasts originates from Bates and Granger [START_REF] Bates | The Combination of Forecasts[END_REF]. Newbold and Granger [START_REF] Newbold | Experience with Forecasting Univariate Time Series and the Combination of Forecasts[END_REF] also suggested combining rules based on variances-covariances of the individual forecasts, while Deutsch et al. [START_REF] Deutsch | The combination of forecasts using changing weights[END_REF] achieved substantially smaller squared forecasts errors combining forecasts with changing weights. Harvey [START_REF] Harvey | Forecasting, structural time series models and the Kalman filter[END_REF] and Hamilton [START_REF] Hamilton | Time series analysis[END_REF] both propose using state space modeling, such as Kalman Filter, for representing dynamic systems where unobserved variables (socalled 'state' variables) can be integrated within an 'observable' model. Finally, Terui and Van Dijk [START_REF] Terui | Combined forecasts from linear and nonlinear time series models[END_REF] also suggest that the combined forecasts perform well, especially with time varying coefficients.

3

US Unemployment dataset

In this application, we forecast the monthly percentage change of the US unemployment rate (UNEMP), as provided by the online Federal Reserve Economic Data (FRED) database of the Federal Reserve Bank of St. Louis 1 . The forecasting perfor-mance of our models is explored over the period of 1972 to 2012, using the last seven years for out-of-sample evaluation. The time series is seasonally adjusted and it is divided into three sub-periods as summarized in Table 1 below: The following graph presents the US unemployment rate for the period under study:

Fig. 1. The US Unemployment Rate

In the absence of any formal theory behind the selection of the inputs of a NN, we conduct some NN experiments and a sensitivity analysis on a pool of potential inputs in the in-sample dataset in order to help our decision. Finally, we select as inputs sets of autoregressive terms of UNEMP that provide the best statistical performance for each network in the test sub-period. These sets are presented in Table 2 below: 50 states or the District of Columbia, who do not reside in institutions (e.g., penal and mental facilities, homes for the aged) and who are not on active duty in the Armed Forces. 1-Jan-75 1-Jan-78 1-Jan-81 1-Jan-84 1-Jan-87 1-Jan-90 1-Jan-93 1-Jan-96 1-Jan-99 1-Jan-02 1-Jan-05 1-Jan-08 In this paper an ARMA model is used to benchmark the efficiency of the NNs' statistical performance. Using as a guide the correlogram and the information criteria in the in-sample subset, we have chosen a restricted ARMA (7, 7) model, where all the coefficients are significant at the 95% confidence interval. The selected ARMA model is presented in equation (1) below:

Y Y Y Y                      (1)
where ˆt Y is the forecasted percentage change of the US unemployment rate.

Neural Networks (NNs)

Neural networks exist in several forms in the literature. The most popular architecture is the Multi-Layer Perceptron (MLP). A standard neural network has at least three layers. The first layer is called the input layer (the number of its nodes corresponds to the number of explanatory variables). The last layer is called the output layer (the number of its nodes corresponds to the number of response variables). An intermediary layer of nodes, the hidden layer, separates the input from the output layer. Its number of nodes defines the amount of complexity the model is capable of fitting. In addition, the input and hidden layer contain an extra node called the bias node. This node has a fixed value of one and has the same function as the intercept in traditional regression models. Normally, each node of one layer has connections to all the other nodes of the next layer.

MLP RNN PSN

UNEMP (1)* UNEMP (1) UNEMP ( 1)

UNEMP (2) UNEMP (3) UNEMP (2)
UNEMP ( 4) UNEMP ( 4) UNEMP ( 3) UNEMP ( 5) UNEMP ( 6) UNEMP ( 6) UNEMP ( 7) UNEMP ( 7) UNEMP ( 8) UNEMP ( 10) UNEMP ( 9) UNEMP ( 10) UNEMP ( 11) UNEMP ( 11) UNEMP ( 11) UNEMP ( 12) UNEMP ( 12) UNEMP ( 12)

The network processes information as follows: the input nodes contain the value of the explanatory variables. Since each node connection represents a weight factor, the information reaches a single hidden layer node as the weighted sum of its inputs. Each node of the hidden layer passes the information through a nonlinear activation function and passes it on to the output layer if the calculated value is above a threshold.

The training of the network (which is the adjustment of its weights in the way that the network maps the input value of the training data to the corresponding output value) starts with randomly chosen weights and proceeds by applying a learning algorithm called backpropagation of errors [START_REF] Shapiro | A Hitchhiker's guide to the techniques of adaptive nonlinear models[END_REF].The learning algorithm simply tries to find those weights which minimize an error function (normally the sum of all squared differences between target and actual values). Since networks with sufficient hidden nodes are able to learn the training data (as well as their outliers and their noise) by heart, it is crucial to stop the training procedure at the right time to prevent overfitting (this is called 'early stopping'). This can be achieved by dividing the dataset into three subsets respectively called the training and test sets used for simulating the data currently available to fit and tune the model and the validation set used for simulating future values. The training of a network is stopped when the mean squared forecasted error is at minimum in the test-sub period. The network parameters are then estimated by fitting the training data using the above mentioned iterative procedure (backpropagation of errors). The iteration length is optimised by maximising the forecasting accuracy for the test dataset. Then the predictive value of the model is evaluated applying it to the validation dataset (out-of-sample dataset).

The Multi-Layer Perceptron Model (MLP)

MLPs are feed-forward layered NN, trained with a back-propagation algorithm. According to Kaastra and Boyd [START_REF] Kaastra | Designing a Neural Network for Forecasting Financial and Economic Time Series[END_REF], they are the most commonly used types of artificial networks in financial time-series forecasting. The training of the MLP network is processed on a three-layered architecture, as described above.

The Recurrent Neural Network (RNN)

The next NN architecture used in this paper is the RNN. For an exact specification of recurrent networks, see Elman [START_REF] Elman | Finding Structure in Time[END_REF]. A simple recurrent network has an activation feedback which embodies short-term memory. The advantages of using recurrent networks over feed-forward networks for modeling non-linear time series have been well documented in the past. However, as mentioned by Tenti [START_REF] Tenti | Forecasting foreign exchange rates using recurrent neural networks[END_REF], "the main disadvantage of RNNs is that they require substantially more connections, and more memory in simulation than standard back-propagation networks" (p. 569), thus resulting in a substantial increase in computational time.

The Psi-Sigma Neural Network (PSN)

The PSNs are a class of Higher Order Neural Networks with a fully connected feedforward structure. Ghosh and Shin [START_REF] Ghosh | The Pi-Sigma Network: An efficient Higher-order Neural Networks for Pattern Classification and Function Approximation[END_REF] were the first to introduce the PSN, trying to reduce the numbers of weights and connections of a Higher Order Neural Network. Their goal was to combine the fast learning property of single-layer networks with the mapping ability of Higher Order Neural Networks and avoid increasing the required number of weights. The price for the flexibility and speed of Psi Sigma networks is that they are not universal approximators. We need to choose a suitable order of approximation (or else the number of hidden units) by considering the estimated function complexity, amount of data and amount of noise present. To overcome this, our code runs simulations for orders two to six and then it presents the best network. The evaluation of the PSN model selected comes in terms of trading performance. 2

5

Forecast Combination Techniques

Kalman Filter

Kalman Filter is an efficient recursive filter that estimates the state of a dynamic system from a series of incomplete and noisy measurements. The time-varying coefficient combination forecast suggested in this paper is shown below:

Measurement Equation:

  3 2 1 , ~0, t t t i i t t i c NNs f a f NID         (2)
State Equation 1 2 , ~(0, )

t t i i t t n a a n n NID    (3)
Where:

 t c
f f f f      (4)
From the above equation we note that the Kalman filtering process also favors the PSN model, which is the model that performs best individually.

Support Vector Regression (SVR)

Vapnik [START_REF] Vapnik | The nature of statistical learning theory[END_REF] established Support Vector Regression (SVR) as a robust technique for constructing data-driven and non-linear empirical regression models. They provide global and unique solutions and do not suffer from multiple local minima (Suykens [START_REF] Suykens | Weighted least squares support vector machines: robustness and sparse approximation[END_REF]). They also present a remarkable ability of balancing model accuracy and model complexity (and Lu et al. [START_REF] Lu | Financial time series forecasting using independent component analysis and support vector regression[END_REF]). The SVR function can be specified as:

( ) ( ) T f x w x b    (5)
where w and b are the regression parameter vectors of the function and φ(x) is the non-linear function that maps the input data vector x into a feature space where the training data exhibit linearity. The ε-sensitive loss Lε function finds the predicted points that lie within the tube created by two slack variables * , i i

 : 0 | ( ) | ( ) , | ( ) | i i i i i if y f x Lx y f x if other               (6)
In other words ε is the degree of model noise insensitivity and Lε finds the predicted values that have at most ε deviations from the actual obtained values yi. The goal is to solve the following argument3 :

Minimize 2 * 1 1 ( ) 2
n i i i C w     subject to * 0 0 0 i i C                and * ( )
() T i i i T i i i y w x b w x b y                           (7)
The above quadratic optimization problem is transformed in a dual problem and its solution is based on the introduction of two Lagrange multipliers * , ii aa and mapping with a kernel function ( ,

) i K x x : * 1 ( ) ( ) ( , ) n i i i f x a a K x x b      where * 0, ii a a C   (8)
Support Vectors (SVs) are called all the xi that contribute to equation ( 8), thus they lie outside the ε-tube, whereas non-SVs lie within the ε-tube. Increasing ε leads to less SVs' selection, whereas decreasing it results to more 'flat' estimates. The norm )

n i i i          
is the training error, as specified by the slack variables. Consequently the introduction of the parameter C satisfies the need to trade model complexity for training error and vice versa (Cherkassky and Ma [START_REF] Cherkassky | Practical selection of SVM parameters and noise estimation for SVM regression[END_REF]). In our application, the NN forecasts are used as inputs for a ε-SVR simulation. A RBF kernel4 is selected and the parameters have been optimized based on cross-validation in the in-sample dataset (ε=0.06, γ= 2.47 and C=0.103).

Statistical Performance

As it is standard in literature, in order to evaluate statistically our forecasts, the RMSE, the MAE, the MAPE and the Theil-U statistics are computed (see Dunis and Williams [START_REF] Dunis | Modelling and Trading the EUR/USD Exchange Rate: Do Neural Network Models Perform Better?[END_REF]). For all four of the error statistics retained the lower the output, the better the forecasting accuracy of the model concerned. In Table 3 we present the statistical performance of all our models in the in-sample period.

Table 3. Summary of the In-Sample Statistical Performance

We note that from our individual forecasts, the PSN statistically outperforms all other models. Both forecast combination techniques improve the forecasting accuracy, but SVR is the superior model regarding all four statistical criteria. Table 4 below summarizes the statistical performance of our models in the out-of-sample period. 4 The RBF kernel equation is From the results above, it is obvious that the statistical performance of the models in the out-of-sample period is consistent with the in-sample one and their ranking remains the same. All NN models outperform the traditional ARMA model. In addition, the PSN outperforms significantly the MLP and RNN in terms of statistical accuracy. The idea of combining NN unemployment forecasts seems indeed very promising, since both Kalman Filter and SVR present improved statistical accuracy also in the out-of-sample period. Moreover, SVR confirms its forecasting superiority over the individual architectures and the Kalman Filter technique. In other words, SVR's adaptive trade-off between model complexity and training error seems more effective than the recursive ability of Kalman Filter to estimate the state of our process.

Concluding Remarks

The motivation for this paper is to investigate the efficiency of a Neural Network (NN) architecture, the Psi Sigma Network (PSN), in forecasting US unemployment and compare the utility of Kalman Filter and Support Vector Regression (SVR) in combining NN forecasts. An Autoregressive Moving Average model (ARMA) and two different NN architectures, a Multi-Layer Perceptron (MLP) and a Recurrent Network (RNN), are used as benchmarks. The statistical performance of our models is estimated throughout the period of 1972-2012, using the last seven years for out-ofsample testing.

As it turns out, the PSN outperforms its benchmark models in terms of statistical accuracy. It is also shown that all the forecast combination approaches outperform in the out-of-sample period all our single models. All NN models beat the traditional ARMA model. In addition, the PSN outperforms significantly the MLP and RNN in terms of statistical accuracy. The idea of combining NN unemployment forecasts seems indeed very promising, since both Kalman Filter and SVR present improved statistical accuracy also in the out-of-sample period. SVR confirms its forecasting superiority over the individual architectures and the Kalman Filter technique. In other, SVR's adaptive trade-off between model complexity and training error seems more effective than the recursive ability of Kalman Filter to estimate the state of our process. The remarkable statistical performance of SVR allows us to conclude that it can be considered as an optimal forecast combination for the models and time series under

term 2 w 1 (

 21 characterizes the complexity (flatness) of the model and the term *

Table 1 .

 1 The US Unemployment Dataset -Neural Networks' Training Dataset

US Unemployment Rate 1 st January 1972 -1 st December 2012

	PERIODS	MONTHS	START DATE	END DATE
	Total Dataset	492	01//01/1972	01/12/2012
	Training Dataset (In-sample)	324	01//01/1972	01/12/1998
	Test Dataset (In-sample)	84	01/01/1999	01/12/2005
	Validation Dataset (Out-of-sample)	84	01/01/2006	01/12/2012

1-Jan-11

Table 2 .

 2 Neural Networks' Inputs

*UNEMPL UNEMP (1) is the first autoregressive term of the UNEMP series 4 Forecasting Models 4.1 Auto-Regressive Moving Average Model (ARMA)

Table 4 .

 4 Summary of the Out-of-sample Statistical Performance

		ARMA	MLP	RNN		PSN	Kalman Filter	SVR
	MAE	1.9941	0.0078	0.0077		0.0073	0.0067	0.0062
	MAPE	65.25%	52.78%	50.17%	47.73%	45.76%	41.52%
	RMSE	2.5903	1.0671	0.9572		0.9045	0.8744	0.8256
	Theil-U	0.6717	0.6142	0.5827		0.5325	0.5017	0.4549
			( , ) exp( K x x 	 	x x 	2	),		0  .
			i			i

  Finally, the results confirm the existing literature, which suggests that nonlinear models, such as NNs, can be used to model macroeconomic series.

	study.
		ARMA	MLP	RNN	PSN	Kalman Filter	SVR
	MAE	0.0332	0.0072	0.0071	0.0061	0.0057	0.0051
	MAPE	67.45%	50.17%	48.97%	44.38%	40.21%	34.33%
	RMSE	2.4043	0.9557	0.9354	0.8927	0.8549	0.8005
	Theil-U	0.5922	0.5654	0.5591	0.4818	0.4657	0.4154

Based on the description given by FRED, the US unemployment rate or civilian unemployment rate represents the number of unemployed as a percentage of the labour force. Labour force data are restricted to people 16 years of age and older, who currently reside in 1 of the

For a complete description of all the neural network models we used and their complete specifications see Sermpinis et al.[START_REF] Sermpinis | Modelling and trading the realised volatility of the FTSE100 futures with higher order neural networks[END_REF].

For a detailed mathematical analysis of the SVR solution see Vapnik[START_REF] Vapnik | The nature of statistical learning theory[END_REF].


Marcio Mendonça
email: mendonca@utfpr.edu.br

Douglas Mastumoto
email: douglas.matsumoto@gmail.com

Lucia V R Arruda
email: lvrarruda@utfpr.edu.br

Elpinik Papageorgiou³
email: epapageorgiou@teilam.gr

Intelligent Systems Applied to the Control of an Industrial Mixer

Keywords: fuzzy cognitive maps, hebbian learning, Arduino platform, process control, fuzzy logic

This paper presents the application of intelligent techniques to control an industrial mixer. Control design is based on hebbian evolution of fuzzy cognitive maps. In this context, this paper develops a dynamical fuzzy cognitive map (D-FCM) based on Hebbian Learning algorithms. Two strategies to update FCM weights are derived. Finally, the D-FCM is used to control an industrial mixer. Simulation results of this control are presented. Additionally, results are provided extending some of the algorithms into the Arduino platform in order to acknowledge the performance of the codes reported in this paper.

Introduction

Artificial Intelligence (AI) has applications in various areas of knowledge, such as mathematical biology, neuroscience, computer science and others. The research area of intelligent computational systems aims to develop methods that try to mimic or approach the capabilities of humans to solve problems. These news methods are looking for emulate human's abilities to cope with very complex processes, based on inaccurate and/or approximated information. However, this information can be obtained from the expert's knowledge and/or operational data or behavior of an industrial system [START_REF] Passino | Fuzzy control[END_REF].

In this context, Fuzzy Cognitive Map (FCM) is a tool for modeling the human knowledge. It can be obtained through linguistic terms, inherent to fuzzy systems, but with a structure similar to the Neural Networks (NN), which facilitates data processing, and has capabilities for training and adaptation. FCM is a technique based on the knowledge that inherits characteristics of Cognitive Maps and Artificial Neural Networks [START_REF] Kosko | Fuzzy cognitive maps[END_REF], [START_REF] Glykas | Fuzzy Cognitive Maps: Advances in Theory, Methodologies, Tools and Applications[END_REF], [START_REF] Kosko | Neural networks and fuzzy systems: a dynamical systems approach to machine intelligence[END_REF], with applications in different areas of knowledge [START_REF] Dickerson | Virtual Worlds as Fuzzy cognitive maps[END_REF], [START_REF] Lee | A cognitive map simulation approach to adjusting the design factors of the electronic commerce web sites[END_REF], [START_REF] Papageorgiou | Novel for supporting medical decision making of different data types based on Fuzzy Cognitive Map Framework[END_REF], [START_REF] Papageorgiou | Combined Fuzzy cognitive map and decision trees model for medical decision making[END_REF], [START_REF] Huang | Application of Fuzzy causal networks to waste water treatment plants[END_REF]. Besides the advantages and characteristics inherited from these primary tech-niques, FCM was originally proposed as a tool to build models or cognitive maps in various fields of knowledge. It makes the tool easier to abstract the information necessary for modeling complex systems, which are similar in the construction to the human reasoning. Thus, FCM aggregates benefits of the acquisition, processing and adaptability from data and information system to be modeled, with a capacity for intelligent decision making due to its heuristic nature. However, FCM has some limitations, especially in time modeling, restricting for applications where the causes and effects occur simultaneously.

In order to circumvent these drawbacks, dynamical fuzzy cognitive maps (D-FCM), can be developed which have the capability to model and manage behaviors of non-linear time-dependent system and often in real time. Examples of different D-FCMs can be found in the recent literature, as examples, we can cite [START_REF] Mendonça | A dynamic fuzzy cognitive map applied to chemical process supervision[END_REF], [START_REF] Miao | Transformation of cognitive maps[END_REF], [START_REF] Papageorgiou | Learning Algorithms for Fuzzy Cognitive Maps[END_REF], [START_REF] Miao | Dynamical cognitive network -an Extension of fuzzy cognitive[END_REF].

Specifically, the work of Mendonça and collaborators [START_REF] Mendonça | A dynamic fuzzy cognitive map applied to chemical process supervision[END_REF] presents a type of D-FCM, which aggregates the occurrence of events and other facilities that makes appropriate this type of cognitive map, for the development of intelligent control and automation in an industrial environment.

In this paper, we use the same D-FCM proposes in [START_REF] Mendonça | A dynamic fuzzy cognitive map applied to chemical process supervision[END_REF] to control an industrial mixing tank. Different from [START_REF] Mendonça | A dynamic fuzzy cognitive map applied to chemical process supervision[END_REF], we use a hebbian algorithm to dynamically adapt the D-FCM weights. In order to validate our D-FCM controller, we compared its performance with a fuzzy logic controller. This comparison is carried out with simulated data. Moreover, to shown the control portability, we embedded the D-FCM controller into a low cost platform based on Arduino.

Development

To demonstrate the evolution of the proposal technique (D-FCM) we will use a case study well known in the literature as seen in [START_REF] Glykas | Fuzzy Cognitive Maps: Advances in Theory, Methodologies, Tools and Applications[END_REF], [START_REF] Axelrod | Structure of decision: the cognitive maps of political elites[END_REF] and [START_REF] Stylios | An Fuzzy Cognitive Maps Approach to Process Control Systems[END_REF] to test level controllers. This case was selected to illustrate the need for refinement of a model based on FCM built exclusively with knowledge. The process shown in Figure 1 consists of a tank with two inlet valves for different liquids, a mixer, an outlet valve for removal of liquid produced by mixing and density meter that measures the quality of the produced liquid. Valves (V 1 ) and (V 2 ) insert two different liquids in the tank. During the reaction of the two liquids, a new liquid characterized by its density value is produced. At this time the valve (V 3 ) empties the tank in accordance with a campaign output flow, but the liquid mixture should be in the specified levels.

Although relatively simple, this process is a TITO (two inputs two outputs) type with coupled variables. To establish the quality of the control system of the produced fluid, a weighting machine placed in the tank measures the (specific gravity) produced liquid.

When the value of the measured variable G (liquid mass) reaches the range of values between the maximum and minimum [Gmin, Gmax] specified, the desired mixed liquid is ready. The removal of liquid is only possible when the volume (V) is in a specified range between the values [Vmin and Vmax]. The control consists to keep these two variables in their operating ranges, as,

V min < V < V max (1)
and

G min < G < G max. (2)
In this study we tried to limit these values from 800 to 850 [mg] for the mass and 850 to 900 [ml] for the volume. According to Papageorgiou et al. [START_REF] Papageorgiou | Fuzzy cognitive maps learning using particle swarm optimization[END_REF], through the observation and analysis of operation of the process is possible for experts to define a list of key concepts related to physical quantities involved. The concepts and cognitive model is based on a known FCM model [START_REF] Papageorgiou | Fuzzy cognitive maps learning using particle swarm optimization[END_REF], having the following concepts and structure: Considering the initial proposed evolution for FCM we will use a D-FCM to control the mixer which should maintain levels of volume and mass within specified limits.


The process model uses the mass conservation principle to derive a set differential equations representing the process used to test the D-FCM controller. As a result the tank volume is the volume over the initial input flow of the intake valves V1 and V2 minus the outflow valve V3. Similarly, the mass of the tank follows the same principle as shown below. The values used for me1 and me2 were 1.0 and 0.9, respectively.

-

) )

The development of the D-FCM is accomplished through two distinct stages. First the D-FCM is developed as a classic FCM where concepts and causal relationships are identified. The concepts can be variables and/or control actions, as already mentioned. However, the heuristic is related to the control condition of volume or weight of the mixture increase, where the inlet valves are closed, thus making it possible to assign inverse causal relationship between the concepts of levels and outlet valves. The output valve defines a positive relationship, when the output flow increases, according to a desired campaign. The flow of the intake valves also increases proportionally. The initial setting is done by using an algorithm based on heuristic optimization method of Simulated Annealing [START_REF] Ghazanfari | Comparing simulated annealing and genetic algorithm in learning fcm[END_REF], in which an initial solution is cast as an initial guess and then solutions with a certain degree of randomness is tested systematically, until the system get the desired response. The figure 2 shows the schematic graph of a D-FCM controller.

Fig. 2. D-FCM Controller

The second stage of development of the D-FCM is responsible for tuning or refinement of the model for dynamic response of the controller. In this case, when a change of output set-point in the campaign occurs, the weights of the causal relationships are tuned. In this paper two options are showed. First, the switching of the causal relationships in different operating points (set-points). It is shown in Table 1 and require different settings offline. One classic genetic algorithm used for tune off-line with 20 individuals, simple crossing and 1% of mutation, however obtained results similar to the method cited. The second option is online and is in accordance with equations 6-11. These equations are based on Hebb learning that increments or decrements the causal relations in accordance with changes in their associated concepts. To perform this function a new kind of concept and relations were included in the cognitive model. The relation selection assigns rules event-driven, which in this case was the change of set-point, levels of volume and weight of the liquid mixture.

The D-FCM uses the concept of selection for switching the set of causal relations according to the basic rules of relation selection. D-FCM works similarly to a DT-FCM (Decision Tree -FCM) [START_REF] Dickerson | Virtual Worlds as Fuzzy cognitive maps[END_REF], in another way, the D-FCM resembles a hybrid tool between a FCM and a temporal cyclic state machine, for example the work of Acampora and Loia [START_REF] Acampora | On the Temporal Granularity in Fuzzy Cognitive Maps[END_REF], and switched by triggering events. The results of the weight of the causal relations by experts through observation of the dynamic behavior of the system by changing the cause and effect of causal relationships are shown in Table 1.

Table 1.

Casual relationship weights In order to establish a correlation and a future comparison between techniques of intelligent systems, a fuzzy controller was also developed. The fuzzy rule base implements weights assignment using the same heuristic control strategy.

As a result, both controllers fuzzy and D-FCM can be run under equal conditions supporting performance comparisons. As an example, we can mention some rules extracted from rule base established by experts:

If VOLUME is HIGH, then V1 is LOW, V2 is LOW; If VOLUME is MEDIAN, then V1 is MEDIAN, V2 is MEDIAN; If WEIGHT is LOW, then V1 is HIGH, V2 is HIGH;
In order to dynamically adapt the D-FCM weights we used the hebbian learning algorithm for FCM that is an adaptation of the classic hebbian method [START_REF] Kosko | Fuzzy cognitive maps[END_REF]. Different proposals and variations of this method applied in tuning or in learning for FCM are known in the literature [START_REF] Papageorgiou | Learning Algorithms for Fuzzy Cognitive Maps[END_REF]. In this paper, the method is used to update the intensity of causal relationships in a deterministic way according to the variation or error in the intensity of the concept or input variable. Specifically, the application of Hebb learning provides control actions as follows: if the weight or volume of the liquid mix increases, the intake valves have a causal relationship negatively intensified and tend to close more quickly. Conversely, if the volume or weight mixture decreases, the intake valves have a causal relationship positively intensified. The mathematical equation is presented in (5).

) ) (5)
∆A i is the concept variation resulting from causal relationship, and it is given by ∆A i = A i (k)-A i (k-1), γ is the learning rate at iteration k.

Causal relationships that have negative causality have negative sign similarly to positive causal relationships. Equations 6-11 show the implementation of the proposal.

)

)

)

)

)

)

Two variations of the Hebbian learning will be presented. In the first case, a global error is assumed as the sum of the variation only of the volume concept. In the second case, the variation of volume, weight and output valve affect the causal relationships (fig. 2), for example, W 51 and W 52 are tuned according to the weight variation ∆A i =V 3 (k)-V 3 (k-1). All values of forgetting factor γ were empirical. Finally, so that variations of the weights had the dynamics needed, any errors or variations of A were multiplied by a factor 20 concerning the gain value.

Experimental Results

The results of D-FCM by Hebbian learning with only the volume parameter variation are shown in Fig. 3, which shows the behavior of the controlled variables within the predetermined range for the volume and the weight of the mixture. It is noteworthy that the controller keeps the variables in the control range and pursues a trajectory according to a campaign where output flow is also predetermined. In this initial experiment, a campaign with a sequence of values ranging from 7, 5 and 10 ml/min can be seen as a set-point output flow (outlet valve). Fig. 4 shows the evolution of the weights of the causal relationships during the process. Fig. 8 shows the results obtained with the Arduino platform providing data of the actuators Valve 1 and Valve 2 with Matlab performing data acquisition. The algorithm switches the sets of causal relations that operate similarly to a DT-FCM (decision tree -FCM), where the activation rules and weights are shown in Table 1.

Similarly, fig. 9 shows the results obtained with the Hebbian learning algorithm for FCM with the three parameters of ΔA i .

Conclusions

Performing a comparison of the results, we see a decrease in the control range in the cases where there is variation in causal relationships by Hebbian learning algorithms. Figures 7 and8 show results.

From the data obtained from Arduino by the variations of the D-FCM embedded in the platform, it is observed that the controlled variables are in well behaved ranges, which suggests that the D-FCM codes have low computational complexity due to the simplicity of mathematical processing compared with the classical fuzzy logic, for example. Thus, we can emphasize the portability and the possibility of developing D-FCM controllers on low cost platforms.

Future studies address a comparison with classical PID controllers, weighted Fuzzy controller or other intelligent controller. Finally, our work will be directed to the application of a real mixer controller into a real data environment.

Fig. 1 .

 1 Fig. 1. Mixer Tank (Source: adapted from Stylios, Groumpos, Georgopoulos, 1999)

Concept 1 -

 1 State of the valve 1 (closed, open or partially open).  Concept 2 -State of the valve 2 (closed, open or partially open).  Concept 3 -State of the valve 3 (closed, open or partially open).  Concept 4 -quantity of fluid (volume) in the tank, which depends on the operational state of the valves V1, V2 and V3.  Concept 5 -value measured by the G sensor for the density of the liquid.

Fig. 3 .

 3 Fig. 3. Inputs and outputs valves, Volume and weight (HL global variation).

Fig. 4 .

 4 Fig. 4. Causal Relationships in the process (HL global variation).

Fig. 5 .Fig. 6 .

 56 Fig. 5. -Inputs and outputs valves, Volume and weight (HL punctual variation)

Fig. 7 .

 7 Fig. 7. Cycle and image of the communication Matlab -Arduino

Fig. 8 .Fig. 9 .

 89 Fig. 8. -Inputs and outputs valves, Volume and weight (DT, Arduino)