Olivier Le Meur
email: olemeur@irisa.fr

Thierry Baccino
email: thierry.baccino@univ-paris8.fr

Methods for comparing scanpaths and saliency maps: strengths and weaknesses

In this paper, we are interested in the computational modelling of visual attention. We report methods commonly used to assess the performance of these kinds of models. We survey the strengths and weaknesses of common assessment methods based on diachronic eye tracking data. We then illustrate the use of some methods to benchmark computational models of visual attention.

Introduction

Eye-tracking is a well-known technique for analysing visual perception and attention shift, and assessing user interfaces. However, until now data analysis from eye-tracking studies has focused on synchronic indicators such as fixation (duration, number, etc) or saccade (amplitude, velocity, etc) rather than diachronic indicators (scanpaths or saliency maps). Synchronic means that an event occurs at a specific point in time, while diachronic means that this event is taken into account over time. We focus in this paper on diachronic measures, and review different ways of analysing sequences of fixations represented as scanpaths or saliency maps.

Visual scanpaths depend on bottom-up and top-down factors such as the task users are asked to perform [START_REF] Simola | Using hidden Markov model to uncover processing states from eye movements in information search tasks[END_REF], the nature of the stimuli [START_REF] Yarbus | Eye movements and vision[END_REF] and the intrinsic variability of subjects [START_REF] Viviani | Eye movements in visual search: cognitive, perceptual and motor control aspects[END_REF]. Being able to measure the difference (or similarity) between two visual behaviours is fundamental both for differentiating the impact of different factors and for understanding what govern our cognitive processes. It also plays a key role in assessing the performance of computational models on overt visual attention, by, for example, evaluating how well saliency-based models predict where observers look.

In this study, we survey common methods for evaluating the difference/similarity between scanpaths and between saliency maps. In the first section, we describe state-of-the-art methods commonly used to compare visual scanpaths. We then consider the comparison methods which involve either two saliency maps or one saliency map plus a set of visual fixations. We first define how human saliency maps are computed, and list some of their most important properties. The strengths and weaknesses of each method are emphasized. In the fourth section, we address inter-observer variability, which reflects the natural dispersion of fixations existing between observers watching the same stimuli. It is important to understand the mechanisms underlying this phenomenon since this dispersion can be used as an upper bound for a prediction. To illustrate the latter point, and the use of common similarity metrics, we compare ground truth and model-predicted saliency maps. Finally, we broaden the scope of this paper by raising a fundamental question: do all visual fixations have the same meaning and role, and is it possible to classify fixations as being bottom-up, cognitive, top-down, or semantic?

Methods for comparing scanpaths

Different metrics are available for comparing two scanpaths, using either distance-based methods (string edit technique or Mannan distance) or vector-based methods. Distance-based methods compare scanpaths only from their spatial characteristics, while vector-based approaches perform the comparison across different dimensions (frequency, time etc). These metrics are more or less complex and relevant depending on the situation to be analysed. However, there is no consensus in the community on the use of a given metric. In this section, we present three metrics: the string edit metric, Mannan's metric and a vector-based metric.

String edit metric

The idea of the string edit metric is that a sequence of fixations on different areas of interest (AOIs) can be translated into a sequence of symbols (numbers or letters) forming strings that are compared. This comparison is carried out by calculating a string edit distance (often called the Levenshtein distance) that gives a measure of the similarity of the strings [START_REF] Levenshtein | Binary codes capable of correcting deletions, insertions and reversals[END_REF]. This technique was originally developed to account for the edit distance between two words, and the measured distance is the number of deletions, insertions or substitutions that are necessary for the two words to be identical (which is also called the alignment procedure). This metric takes as input two strings (coding AOIs) and computes the minimum number of edits needed to transform one string into the other. A cost is associated with each transformation and each character. The goal is to find a series of transformations that minimizes the cost of aligning one sequence with the other. The total cost is the edit distance between the two strings. When the cost is minimal, the similarity between the two strings is maximal (i.e. when the two strings are identical, the distance is equal to 0). Conversely, the distance increases with the cost and therefore with the dissimilarity between the two strings. Figure 1 illustrates the method. The Levenshtein distance is the most common way to compare scanpaths [START_REF] Josephson | Attention to repeated images on the World-Wide Web: Another look at scanpath theory[END_REF][START_REF] Privitera | Algorithms for defining visual regions-of-interest: Comparison with eye fixations[END_REF] and has been widely used for assessing the usability of web pages [START_REF] Baccino | La Lecture électronique [Digital Reading[END_REF]. The string edit distance can be computed using a dynamic programming technique (the WagFish algorithm [START_REF] Wagner | The string-to-string correction problem[END_REF]) that incrementally computes optimal alignments (minimizing the cost). The Levenshtein distance is not the only string edit distance that can be used for scanpaths. Others are described below:

• LCS is the length of the Longest Common Subsequence, which represents the score obtained by allowing only addition and deletion, not substitution; • Damerau-Levenshtein distance allows addition, deletion, substitution and the transposition of two adjacent characters; and • Hamming distance only allows substitution (and hence, only applies to strings of the same length).

The advantage of the string edit technique is that it is easily computed and keeps the order of fixations. It is also possible to compare observed scanpaths to predicted scanpaths when certain visual profiles are expected from the cognitive model used by the researcher [START_REF] Chanceaux | Towards a model of information seeking by integrating visual, semantic and memory maps[END_REF]. However, several drawbacks have to be underlined:

• Since the string edit is based on a comparison of the sequence of fixations occurring in pre-defined AOIs, the question is how to define these AOIs. There are two ways: automatically gridded AOIs or content-based AOIs. The former is built by putting a grid of equally sized areas across the visual material, but for the latter the meaningful regions of the stimulus need to be subjectively chosen. Whatever AOIs are constructed, the string edit method means that only the quantized spatial position of the visual fixations are taken into account. Hence, some small differences in scanpaths may change the string while others produce the same string. • The string edit method is limited when certain AOIs have not been fixated so there is a good deal of missing data.

Mannan's metric

The Mannan distance [START_REF] Mannan | Automatic control of saccadic eye movements made in visual inspection of briefly presented 2-D images[END_REF], 1996[START_REF] Mannan | Fixation sequences made during visual examination of briefly presented 2D images[END_REF] is another metric comparing scanpaths but based on their spatial properties rather than their temporal dimensions, in the sense that the order of fixations is completely ignored. The Mannan distance compares the similarity between scanpaths by calculating the distance between each fixation in one scanpath and its nearest neighbour in the other scanpath. A similarity index (Is) represents the average linear distance between two scanpaths (D), with randomized scanpaths having the same size (Dr). These randomly generated scanpaths are used for weighting the sequence of real fixations, taking into account the fact that real scanpaths may convey a randomized component. The similarity index (Is) is given by

𝐼𝑠 = �1 - 𝐷 𝐷𝑟 � × 100 D is a measure of distance given by 𝐷² = 𝑛 1 ∑ 𝑑 2𝑗 2 𝑛 2 𝑗=1 2𝑛 1 𝑛 2 (𝑎 2 + 𝑏 2 ) + 𝑛 2 ∑ 𝑑 1𝑖 2 𝑛 1 𝑖=1 2𝑛 1 𝑛 2 (𝑎 2 + 𝑏 2 )
where,

• 𝑛 1 and 𝑛 2 are the number of fixations in the two traces; • 𝑑 1𝑖 is the distance between the i th fixation in the first trace and its nearest neighbour in the second trace; • 𝑑 2𝑗 the distance between the j th fixation in the second trace and that of its nearest neighbour in the second one; • 𝑎 and 𝑏 are the side lengths of the image; and • Dr is the distance between two sets of random locations.

The values returned by the algorithm (𝐼𝑠) range from 0 (random scanpath) to 100 (identity). The main drawbacks of this technique are:

• The Mannan distance does not take into account the temporal order of fixation sequence. This means that two sequences of fixation having a reversed order but with an identical spatial configuration give the same Mannan distance. • A difficult problem occurs when the two scanpaths have very different size (the number of fixations between them is very different). The Mannan distance may show a great similarity while the shapes of the scanpaths are definitely different. The Mannan distance is not tolerant to high variability between scanpaths.

Vector-based metric

An interesting method was recently proposed by [START_REF] Jarodzka | A vector-based, multidimensional scanpath similarity measure[END_REF]. Each scanpath is viewed as a sequence of geometric vectors that corresponds to subsequent saccades of the scanpath. The vector representation shows the length and the direction of each saccade. A saccade is defined by a starting position (fixation n) and ending position (fixation n+1). Then a scanpath with n fixations is represented by a set of n-1 vectors, and several properties can therefore be preserved, such as the shape of the scanpath, the scanpath length, and the position and duration of fixations. The sequences that have to be compared are aligned according to their shapes (although this alignment can be performed on other dimensions: length, durations, angle, etc).

Each vector of one scanpath corresponds to one or more vectors of another scanpath, such that the path in the matrix of similarity between the vectors going from (1, 1) (similarity between the first vectors) to (n,m) (similarity between the last vectors) is the shortest one. Once the scanpaths are aligned, various measures of similarity between vectors (or sequences of vectors) can be used, such as average difference in amplitude, average distance between fixations and average difference in duration.

Figure 2: Alignment using saccadic vectors. The alignment procedure attempts to match the 5 vectors of the two scanpaths. The best match is the following: 1-2/1; 3/2; 4/3; 5/4-5.

For example, Figure 2 shows two scanpaths A and B (the first saccade is going upward). The alignment procedure attempts to match the five vectors (for the five consecutive saccades) of the participant scanpath with the four vectors of the model scanpath. Saccades 1 and 2 of scanpath A are aligned with saccade 1 of scanpath B, saccade 3A is aligned with saccade 2B, etc. Once the scanpaths are aligned, similarity measures are computed for each alignment. Jarodzka's procedure ends up with five measures of similarity (difference in shape, amplitude and direction between saccade vectors, distance between fixations and fixation durations).

This vector-based alignment procedure has a number of advantages over the string edit method. The first is that it does not need to determine pre-defined AOIs (and is therefore not dependent on a quantization of space). The second one is that it can align scanpaths not only on spatial dimension but also on any dimension available in saccade vectors (angle, duration, length, etc). For example, Lemaire, Guérin-Dugué, [START_REF] Lemaire | A cognitive computational model of eye movements investigating visual strategies on textual material[END_REF] used the spatial distance between saccades, the angle between saccades, and the difference of amplitude to realize the alignment. Thirdly, this alignment technique provides more detailed information on the type of (dis)similarity of two scanpaths according to the dimensions chosen. Lastly, the new measure deals with temporal issues, not only fixation durations, but it also successfully deals with shifts in time and variable scanpath lengths. The major drawbacks are:

• This measure only compares two scanpaths. Sometimes the overall aim is to compare whole groups of participants with each other. • It is presumed that fixations and saccades must occur. Other eye movements such as smooth pursuit are not handled. Smooth pursuit movements are important when watching a video. However, the problem may be solved if it is possible to represent smooth pursuit as a series of short vectors which are not clustered into one long vector. • This alignment procedure is independent of the stimulus content. However, the chosen dimensions may be weighted by some semantic values carefully selected by the researcher.

Methods for comparing saliency maps

Comparing two scanpaths requires taking a number of factors, such as the temporal dimension or the alignment procedure, into account. To overcome these problems, another kind of method can be used. In this section, we focus on approaches involving two bi-dimensional maps. We first briefly describe how the visual fixations are used to compute a continuous saliency map. Secondly, we describe three common methods used to evaluate the degree of similarity between two saliency maps: a correlation-based measure, the Kullback-Leibler divergence and ROC analysis.

From a discrete fixation map to a continuous saliency map

A discrete fixation map 𝑓 𝑖 for the i th where 𝑥 is a vector representing the spatial coordinates and 𝑥 𝑓 (𝑘) is the spatial coordinates of the observer is classically defined as

𝑓 𝑖 (𝑥) = � 𝛿(𝑥 -𝑥 𝑓 (𝑘) ) 𝑀 𝑘=1
k th visual fixation. The value M is the number of visual fixations for the For N observers, the final fixation map f is given by:

𝑓(𝑥) = 1 𝑁 � 𝑓 𝑖 (𝑥) 𝑁 𝑖=1 . i th observer. 𝛿(. ) is the Kronecker symbol (𝛿(𝑡) = 1 if t=1, otherwise 𝛿(𝑡) = 0).
A saliency map S is then deduced by convolving the fixation map 𝑓 by an isotropic bi-dimensional Gaussian function

𝑆(𝑥) = 𝑓(𝑥) * 𝐺 𝜎 (𝑥)
where σ is the standard deviation of the Gaussian function. It is commonly accepted that σ should be set to one degree of visual angle. One degree of visual angle represents an estimate of the size of the fovea. The standard deviation depends on the experimental setup (size of the screen and viewing distance). It is also implicitly assumed that a fixation can be approximated by a Gaussian distribution [START_REF] Bruce | Saliency based on information maximisation[END_REF][START_REF] Velichkovsky | Attention and communication: eye-movement-based research paradigms[END_REF]. An example of fixation and saliency maps is given by Figure 3. A heat map which is a simple coloured representation of the continuous saliency map is also shown. Red areas pertain to salient areas whereas blue areas are for non-salient areas. Note that the fixation map illustrated by Figure 3 is not exactly the one defined by the previous formula. Throughout this section, we use the two continuous saliency maps shown in Figure 4 to illustrate the comparison methods. Both maps are obtained from visual fixations of three observers.

Fixation duration is not taken into account in the computation of the continuous saliency map. In 2005, Itti showed that there was no significant correlation between model-predicted saliency and duration of fixation. Fixation duration is often considered to reflect the depth and the speed of visual processing in the brain. The longer the fixation duration, the deeper the visual processing [START_REF] Henderson | Regarding scenes[END_REF][START_REF] Just | Eye fixations and cognitive processes[END_REF]. Total fixation time (the cumulative duration of fixations within a region) can be used to gauge the amount of total cognitive processing engaged with the fixated information [START_REF] Rayner | Eye movements in reading and information processing: 20 years of research[END_REF]. There are a number of factors that influence the duration of fixations. Among these factors, the visual quality of the displayed stimuli plays an important role -as suggested by Mannan et al.'s experiment [START_REF] Mannan | Automatic control of saccadic eye movements made in visual inspection of briefly presented 2-D images[END_REF]. They presented filtered and unfiltered photos to observers and reported a significant increase in the fixation duration for the filtered photos. Another factor is related to the number of objects present in the scene. [START_REF] Irwin | Eye movements and scene perception: memory for things observed[END_REF] observed that the duration of fixations increases with the number of objects in the scene.

Heat map

Saliency map Heat map Saliency map

(a) (b) (c) (d)
Figure 4: Heat maps and continuous saliency maps obtained from fixations of two groups of three observers.

Correlation-based measures

The Pearson correlation coefficient r between two maps 𝐻 and 𝑃 is defined as:

𝑟 𝐻,𝑃 = 𝑐𝑜𝑣(𝐻, 𝑃) 𝜎 𝐻 𝜎 𝑃
where 𝑐𝑜𝑣(𝐻, 𝑃) is the covariance between 𝐻 and 𝑃, and 𝜎 𝐻 and 𝜎 𝑃 represent the standard deviation of maps 𝐻 and 𝑃, respectively.

The linear correlation coefficient has a value between -1 and 1. A value of 0 indicates that there is no linear correlation between the two maps. Values close to zero indicate a poor correlation between the two sets. A value of 1 indicates a perfect correlation. The sign of r is helpful in determining whether data share the same structure. A value of -1 also indicates a perfect correlation, but the data vary together in opposite directions.

This indicator is very simple to compute and is invariant to linear transformation. Several studies have used this metric to assess the performance of computational models of visual attention [START_REF] Jost | Assessing the contribution of color in visual attention[END_REF][START_REF] Le Meur | A coherent computational approach to model bottom-up visual attention[END_REF][START_REF] Rajashekar | Gaffe: A gaze-attentive fixation finding engine[END_REF]. Correlations are usually reported with degrees of freedom (the total population minus 2) in parentheses and the significance level. For instance, the two continuous saliency maps illustrated by Figure 4 are strongly correlated, r(393214)=0.826, p<.001.

Note that the Spearman's rank correlation can also be used to measure the similarity between two sets of data [START_REF] Toet | Computational versus Psychophysical Bottom-Up Image Saliency: A Comparative Evaluation Study[END_REF].

The Kullback-Leibler divergence

The Kullback-Leibler divergence is used to estimate the overall dissimilarity between two probability density functions. Let us define two discrete distributions 𝑅 and 𝑃 with probability density functions 𝑟 𝑘 and 𝑝 𝑘 . The KL-divergence between 𝑅 and 𝑃 is given by the relative entropy of 𝑃 with respect to R:

𝐾𝐿(𝑅, 𝑃) = ∑ 𝑝 𝑘 𝑙𝑜𝑔 𝑟 𝑘 𝑝 𝑘 𝑘 .
The KL-divergence is only defined if 𝑟 𝑘 and 𝑝 𝑘 both sum to 1 and if 𝑟 𝑘 > 0 for any 𝑘 such that 𝑝 𝑘 > 0.

The KL-divergence is not a distance, since it is not symmetric and does not satisfy the triangle inequality. The KL-divergence is non-linear. It varies in the range of zero to infinity. A zero value indicates that the two probability density functions are strictly equal. The fact that the KL-divergence does not have a well-defined upper bound is a strong drawback.

In our context we have to compare two bi-dimensional saliency maps (𝐻 and 𝑃). We first transform these maps into two bi-dimensional probability density functions by dividing each location of the map by the sum of all pixel values. The probability that an observer focuses on position 𝑥 is given by:

𝑝 𝐻 (𝑥) = 𝐻(𝑥) + 𝜖 ∑ (𝐻(𝑖) + 𝜖) 𝑖 𝑝 𝑃 (𝑥) = 𝑃(𝑥) + 𝜖 ∑ (𝑃(𝑖) + 𝜖) 𝑖
where 𝜖 is a small constant to avoid division by zero. If we consider the example of Figure 4, we compute the KL-divergence by first considering the saliency map (b) as the reference and secondly the saliency map (d) as the reference. We obtain KL=3.33 and KL=7.06 respectively. As the KL-divergence is not a distance, the results are not the same. They indicate that the overall dissimilarity is highest when the continuous saliency map (d) is taken as the reference.

Receiver Operating Characteristic Analysis

The Receiver Operating Characteristic (ROC) analysis [START_REF] Green | Signal detection theory and psychophysics[END_REF] is probably the most popular and most widely used method in the community for assessing the degree of similarity of two saliency maps. ROC analysis classically involves two sets of data: the first is from the ground truth (also called the actual values) and the second is the prediction (also called the outcomes).

Here we perform ROC analysis between two maps. It is also common to encounter a second method in the literature that involves fixation points and a saliency map. This method is described in section 4.

Continuous saliency maps are processed as a binary classifier applied on every pixel. It means that the image pixels of the ground truth as well as those of the prediction are classified as fixated (or salient) or as not fixated (or not salient). A simple threshold operation is used for this purpose. However, two different processes are used depending on whether the ground truth or the prediction is considered:

• Thresholding the ground truth Figure 5 : the continuous saliency map is thresholded with a constant threshold in order to keep a given percentage of image pixels. For instance, we can keep the top 2, 5, 10, or 20 % salient pixels of the map, as illustrated by . This threshold is called T x,G (G for the ground truth and x indicating the percentage of image considered as being fixated).

• Thresholding the prediction: the threshold is systematically moved between the minimum and the maximum values of the map. A high threshold value corresponds to an over-detection whereas a smaller threshold affects the most salient areas of the map. This threshold is called T x,P (P for the prediction and x indicating the For each pair of thresholds, four numbers featuring the quality of the classification are computed. They represent the true positives (TP), the false positives (FP), the false negatives (FN) and the true negatives (TN). The true positive number is the number of fixated pixels in the ground truth that are also labelled as fixated in the prediction. 7. A more efficient approximation can be obtained by a trapezoid approximation: rather than computing the area of rectangles, the AUC is given by summing the area of trapezoids. In our example, the AUC value is 0.83.

20% of ground truth

Prediction Classification

Hybrid methods

So far we have focused on similarity metrics involving two scanpaths or two saliency maps. In this section we describe methods based on a saliency map and a set of fixation points. We call this kind of method hybrid as it mixes two types of information. Four approaches are presented: ROC analysis, Normalized Scanpath Saliency, percentile and the Kullback-Leibler divergence.

Receiver Operating Characteristic Analysis

The ROC analysis is performed here between a continuous saliency map and a set of fixations. The method tests how the saliency at the points of human fixation compares to the saliency at non-fixated points. As in the previous section, the continuous saliency map is thresholded to keep a given percentage of pixels of the map. Each pixel is then labelled as either fixated or not-fixated. For each threshold the observer's fixations are laid down on the thresholded map. The true positive (fixations that fall on fixated areas) and the false negative (fixations that fall on non-fixated areas) are determined (as illustrated by Figure 8). A curve that shows the TPR (or hit rate) as a function of the threshold can be plotted. Note that the percentage of the image considered to be salient is in the range of 0 to 100%. If the fixated and non-fixated locations cannot be discriminated, the AUC will be 0.5. This first analysis method is used in papers such as [START_REF] Tatler | Visual correlates of fixation selection: effects of scale and time[END_REF][START_REF] Tatler | Visual correlates of fixation selection: effects of scale and time[END_REF][START_REF] Torralba | Contextual guidance of eye movements and attention in real-world scenes: the role of global features in object search[END_REF][START_REF] Torralba | Contextual guidance of eye movements and attention in real-world scenes: the role of global features in object search[END_REF]. Although interesting, this method is not sensitive to the false alarm rate. To deal with the previous limitation, a set of control points, corresponding to non-fixated points, can be generated. Two methods commonly encountered in the literature are discussed. The first method is the simplest one and consists in selecting control points from either a uniform or a random distribution. This solution does not take into account the fact that fixations are distributed neither evenly nor randomly throughout the scene, as illustrated by Figure 9. The second method, proposed by [START_REF] Einhauser | Does luminance-contrast contribute to a saliency for overt visual attention[END_REF][START_REF] Einhauser | Does luminance-contrast contribute to a saliency for overt visual attention[END_REF][START_REF] Tatler | Visual correlates of fixation selection: effects of scale and time[END_REF], defines control points by choosing locations randomly from a distribution of all fixation locations that occurred at the same time, but on other images. This way to define the control point is important for different reasons. First, as the fixations come from the same observer, so the same bias, systematic tendency or peculiarity of the observer are taken into account. This is illustrated by Figure 10. These factors then have a limited influence on the classification results. Among them, the most important influence is the central bias illustrated in Figure 9. A number of factors can explain this central tendency. The centre might reflect an advantageous viewing position for extracting visual information [START_REF] Renninger | Where to look next? Eye movements reduce local uncertainty[END_REF][START_REF] Tatler | The central fixation bias in scene viewing: selecting an optimal viewing position independently of motor biases and image feature distributions[END_REF]. However, [START_REF] Tatler | The central fixation bias in scene viewing: selecting an optimal viewing position independently of motor biases and image feature distributions[END_REF] noticed that the distribution of low-level visual features over the picture has no significant impact on this bias. In addition, this tendency to look at the centre of images is particularly difficult to remove. Different strategies were tried by [START_REF] Bindemann | Scene and screen center bias early eye movements in scene viewing[END_REF] to reduce this bias, but without success. This author concluded that the screen-based central fixation bias might be an inescapable feature of scene viewing under laboratory conditions. Secondly, the set of control points has to stem from the same time interval as the set that is analysed. Indeed, bottom-up and top-down influences are not similar over time. For instance, bottom-up influences are maximal just after the stimulus onset. Top-down influences tend to increase with viewing time, leading to a stronger dispersion between observers (see section 5). Although the second method is more robust than the first one, the method has a serious flaw. It underestimates the salience of areas which are more or less centred in the image. In a similar fashion to the method in section 3.4, the control points and the fixation points are then used to plot an ROC curve. For each threshold the observer's fixations and the control ones are laid down on the thresholded map. The true positive rate (fixations that fall on fixated areas) and the false positive rate are determined. From this ROC curve the AUC is computed. The confidence interval is computed by using a non-parametric bootstrap technique [START_REF] Efron | An Introduction to the Bootstrap[END_REF]. Many samples having the same size as the original set of human fixations are generated by sampling with replacement. These samples are called bootstrap samples. In general 1,000 bootstrap samples are created. Each bootstrap sample is used as a set of control fixations. The ROC area between the continuous saliency map and the points of human fixation plus the control points is computed. The bootstrap distribution of each ROC analysis is computed and a bootstrap percentile confidence interval is determined by percentiles of the bootstrap distribution, leaving off Sometimes, the quality of the classification relies on the equal error rate (EER). The equal error rate is the location on an ROC curve where the false positive rate and the true positive rate are equal (i.e. the error at which false alarms equal the miss rate FPR=1-TPR). As with the AUC, the EER is used to compare the accuracy of the prediction. In general, the system with the lowest EER is the most accurate.

Normalized Scanpath Saliency

The Normalized Scanpath Saliency (NSS) [START_REF] Peters | Components of bottom-up gaze allocation in natural images[END_REF]) is a metric that involves a saliency map and a set of fixations. The idea is to measure the saliency values at fixation locations along a subject's scanpath.

The first thing to do is to standardize the saliency values in order to have a zero mean and unit standard deviation. It is simply given by 𝑍 𝑆𝑀 (𝑥) =

𝑆𝑀(𝑥) -𝜇 𝜎

where 𝑍 𝑆𝑀 is the standardized saliency map and

𝜇 = 1 |𝐼| � 𝑆𝑀(𝑥 𝑡 ) 𝑡∈𝐼 𝜎 = � 1 |𝐼| �(𝑆𝑀(𝑥 𝑡 ) 𝑡∈𝐼 -𝜇) 2
where the operator | |

. indicates the number of pixels of the picture. For a given coordinate, the quantity 𝑍 𝑆𝑀 (𝑥 𝑖 ) represents the distance between the saliency value at 𝑥 𝑖 and the average of saliency expressed in units of the standard deviation. This value is negative when the saliency value at the fixation locations is below the mean, positive when above. To take account of the fact that we do not focus accurately on a particular point, the NSS value for a given fixation location is computed on a small neighbourhood centred on that location.

𝑁𝑆𝑆 �𝑥 𝑓

(𝑘) � = � 𝐾 ℎ �𝑥 𝑓 (𝑘) -𝑥 𝑗 � 𝑍 𝑆𝑀 (𝑥 𝑗 ) 𝑗∈𝜋
where K is a kernel with a bandwidth ℎ and 𝜋 is a neighbourhood.

The NSS is the average of 𝑁𝑆𝑆 �𝑥 𝑓 (𝑘) � for all fixations 𝑀 of an observer. It is given by

𝑁𝑆𝑆 = 1 𝑀 ∑ 𝑁𝑆𝑆 �𝑥 𝑓 (𝑘) � 𝑀 𝑘=1
.

Figure 11 illustrates the computation of the NSS value for a scanpath composed of 8 visual fixations. In this example, the average NSS value is 0.3, indicating a good correspondence between the model-predicted saliency map and the observer's scanpath.

Percentile

In 2008, Peters and Itti designed a metric called percentile [START_REF] Peters | Applying computational tools to predict gaze direction in interactive visual environments[END_REF]. A percentile value 𝑃 �𝑥 𝑓 (𝑘) � is computed for each location of fixation points 𝑥 𝑓 (𝑘) . This score is the ratio between the number of locations in the saliency map with values smaller than the saliency value at point 𝑥 𝑓 (𝑘) and the set of all locations. The percentile value is defined as follows:

𝑃 �𝑥 𝑓 (𝑘) � = 100 × ��𝑥 ∈ 𝑋: 𝑆𝑀(𝑥) < 𝑆𝑀(𝑥 𝑓 (𝑘) )��

|𝑆𝑀|

where 𝑋 is the set of locations of the saliency map SM and 𝑥 𝑓 (𝑘) is the location of the k th fixation.

The final score is the average of 𝑃 �𝑥 𝑓 (𝑘) � for all fixations of an observer. By definition, the percentile metric has a well-defined upper bound (100%) indicating the highest similarity between fixation points and saliency map. The chance level is 50%.

| |

. indicates set size.

The Kullback-Leibler divergence

The KL-divergence, defined in section 3.3, is used here to compute the dissimilarity between the histogram of saliency sampled at eye fixations and that sampled at random locations. [START_REF] Itti | Bayesian surprise attracts human attention[END_REF] were the first to use this method. The set of control points (or the set of non-fixated points) are drawn from a uniform spatial distribution. However, human fixations are not randomly distributed, since they are governed by various factors such as the central bias explained earlier. To be more agnostic to this kind of mechanism, Zhang, Tong, Marks, Shan & Cottrell in 2008 measured the KL-divergence between the saliency distribution of fixated points of a test image and the saliency distribution at the same pixel locations but of a randomly chosen image from the test set. To evaluate the variability of the score, the evaluation was repeated 100 times with 100 different sets of control points.

Contrary to the previous KL-divergence method of section 3.3, a good prediction has a high KL-divergence score. Indeed, as the reference distribution represents chance, the saliency computed at human-fixated locations should be higher than that computed at random locations.

Measuring a realistic upper-bound

Most of the methods mentioned above have a well-defined theoretical upper-bound. When assessing the performance of a computational model, it is then reasonable to seek to approach this upper-bound. For instance, according to the ROC analysis, an AUC equal or close to one would indicate a very good performance. In our context this goal is almost impossible to reach. Indeed, there is a natural dispersion of fixations among different subjects looking at the same image. This dispersion (also called Inter-Observer Congruency -IOC 1 1 Note that a high dispersion corresponds to a lack of congruency.

) is contingent upon a number of factors. Firstly, [START_REF] Tatler | Visual correlates of fixation selection: effects of scale and time[END_REF] showed that the consistency between visual fixations of different participants is high just after the stimulus onset but progressively decreases over time. Among the reasons that might explain this variability, the most probable one concerns the time course of bottom-up and top-down mechanisms. Just after the stimulus onset our attention is mostly steered by low-level visual features, whereas top-down mechanisms become more influential after several seconds of viewing. The second factor concerns the visual content itself. In the case where there is nothing that stands out from the background, the IOC would be small. On the contrary, a visual scene composed of salient areas would presumably attract our visual attention, leading to high congruency. The presence of particular features such as human faces, human beings or animals tends to increase the consistency between observers. A number of studies have shown that we are able to identify and recognize human faces and animals very quickly in a natural scene [START_REF] Delorme | Key visual features for rapid categorization of animals in natural scenes[END_REF][START_REF] Rousselet | Is it an animal? Is it a human face? Fast processing in upright and inverted natural scenes[END_REF]. Whereas human faces and animals have the ability to attract our attention, decreasing the dispersion between observers, this ability is modulated by novelty or even emotion. Althoff and Cohen's 1999 study is a good example of this point. They investigated the effect of memory or prior experience on eye movements. They found that visual scanpaths made when viewing famous faces were more variable than those made when viewing non-famous faces. A third factor that could account for the variance between people might be related to cultural differences. [START_REF] Nisbett | The geography of thought: how Asians and Westerners think differently[END_REF] compared the visual scan pattern of American and Asian populations. He found that Asian people tend to look more at the background and spend less time on focal objects than American people. However, a recent study casts doubt on the influence of cultural differences on oculomotor behaviour [START_REF] Rayner | Eye movements when looking at unusual-weird scenes: are there cultural differences?[END_REF].

The inter-observer congruency can be measured by using a one-against-all approach, also called "leave one out" [START_REF] Torralba | Contextual guidance of eye movements and attention in real-world scenes: the role of global features in object search[END_REF]. It consists in computing the degree of similarity between the fixations of one observer and those of the other subjects. The final value is obtained by averaging the degree of similarity over all subjects.

In this paper, the ROC metric is used to compute the degree of similarity, as proposed by [START_REF] Torralba | Contextual guidance of eye movements and attention in real-world scenes: the role of global features in object search[END_REF]. The first step consists of building a saliency map from the visual fixations of all observers except one (the i th observer). This map is thresholded so that the most fixated areas are set to 1 and the other areas are set to 0. To assess the degree of similarity between the 4.1

i th observer and the other subjects, the hit rate (as described in section ), i.e. the percentage of fixations that fall into the fixated regions, is computed. Iterating over all subjects and averaging the scores gives the Inter-Observer Congruency. Figure 12 illustrates this method.

Figure 12: IOC measurement. On the left, spatial coordinates of visual fixations for each observer are given. By considering all fixations except those from the i th observer, a heat map is computed (on the right). After an adaptive binarization, we count the number of fixations of the i th observer that fall into salient regions (white regions on the bottom image). The second picture (Stream) is a mountain landscape where no element strongly attracts our attention. For both pictures, the congruency decreases over time as expected. The highest congruency, observed at the beginning of the viewing, is likely to be due to the bottom-up influences and the central bias. It is interesting to emphasize that the congruency for the Parrots image is significantly higher than that observed for the Stream image. This observation is mainly due to the attractiveness of the two contents. When there is nothing in the scene that catches our attention, observers are not 'unconsciously-constrained' and they just do not explore the visual scene in the same way. Assessing the IOC is fundamental to evaluating the performance of the saliency algorithm, although this is overlooked most of the time. An absolute score between a set of fixations and a predicted map is interesting but is not sufficient to draw any conclusions. A low score of prediction does not systematically indicate that the saliency model performs poorly. Such a statement would be true if the dispersion between the observers is low, but false otherwise. Therefore it is much more relevant to compare the performance of computational models to the IOC [START_REF] Judd | Learning to predict where humans look[END_REF][START_REF] Torralba | Contextual guidance of eye movements and attention in real-world scenes: the role of global features in object search[END_REF] or to express the performance directly by normalizing the similarity score by the inter-observer congruency [START_REF] Zhao | Learning a saliency map using fixated locations in natural scenes[END_REF]. The normalized score would be close to 1 for a good prediction.

Example: performance of state-of-the-art computational models

In this section, we examine the performance of the most prominent saliency models that have been proposed in the literature. The quality of the predicted saliency maps is given here by two metrics: the hit rate (see section 4.1) and the NSS (see section 4.2). These metrics are hybrid metrics, since they involve a set of visual fixations and a map. We believe that these metrics are the best way to assess the relevance of a predicted saliency map. Compared to saliency map-based methods, hybrid methods are non-parametric. Human saliency maps are obtained by convolving a fixation map by a 2D Gaussian function, which is parameterized by its mean and its standard deviation. Note that instead of using the hit rate, we could have used an ROC analysis.

To perform the analysis, we use two eye-tracking datasets that are available on the Internet. They are described in section 6.1. We present and comment on each model's performance in section 6.2.

Eye-tracking datasets

Eye-tracking is nowadays a common solution for studying visual perception. Since 2000, some eye-tracking data can be freely downloaded from the Internet for scientific purposes. Table 1 gives the main characteristics of the most important data collections on the web. They are composed of stimuli that represent landscape, outdoor or indoor scenes. Some of them are composed of high-level information such as people, faces, animals and text.

These datasets can be used to evaluate the performance of computational models. There is only an implicit consensus on how to set up an eye-tracking test. There is no document that accurately describes what must be done and what must be avoided in the experimental setting. For instance, should observers perform a task when viewing stimuli or not? Do the methods used to identify fixations and saccades from the raw eye data give similar results? We have to be aware that these datasets have been collected in different environments and with different apparatus and settings. To evaluate the performance of saliency models, it is highly recommended that more than one dataset is used in order to strengthen the findings.

Benchmark

We compare the performance of four state-of-the-art models: Itti's model [START_REF] Itti | A model for saliency-based visual attention for rapid scene analysis[END_REF], Le Meur's model [START_REF] Le Meur | A coherent computational approach to model bottom-up visual attention[END_REF]), Bruce's model [START_REF] Bruce | Saliency, attention, and visual search: An information theoretic approach[END_REF]) and Judd's model [START_REF] Judd | Learning to predict where humans look[END_REF]. (For a brief review of saliency models, see Le [START_REF] Le Meur | What we see is most likely to be what matters: Visual attention and applications[END_REF].

Two eye-tracking datasets (Le Meur and Bruce, see Table 1) are used. The degree of similarity between ground truth and model-predicted saliency is evaluated by using the ROC analysis (hit rate) and the NSS metric. Figure 14 gives the ROC curve indicating the performance of different saliency models averaged over all testing images. The method used here is the method described at the beginning of section 4.1. The upper-bound, i.e. the inter-observer variability, was computed by the method proposed by [START_REF] Torralba | Contextual guidance of eye movements and attention in real-world scenes: the role of global features in object search[END_REF] and described in section 5. Table 2 gives the average NSS value over the two tested datasets.

Under the ROC metric, Judd's model has the highest performance, as illustrated by Figure 14. This result was expected, since this model uses specific detectors (face detection for instance) that improve the ability to detect salient areas. In addition, this model uses a function to favour the centre of the picture in order to take the central bias into account. However, the results are more contrasted under the NSS metric shown in Table 2. On average across both databases, Judd's model is still the highest performing. On Bruce's dataset, Itti's model performs the best, with a value of 0.99, whereas Judd's model performs at 0.87. The model ranking is therefore dependent on the metric used. It is therefore fundamental to use more than one metric when assessing the performance of computational models of visual attention. Current computational models of visual attention focus on identifying fixated locations of salient areas. From an input picture, a model computes a topographic map indicating the most visually interesting parts. This prediction is then compared to ground truth fixations. The evaluation methodology seems to be appropriate. Unfortunately an important point is overlooked. By doing this kind of comparison, most researchers have implicitly assumed that fixations, whatever their durations, saccade amplitudes and start-time, are all similar. In this section, we emphasize the fact that different populations of fixations may exist.

Fixations differ in both their durations and their saccade amplitudes during real-world scene viewing. Antes in 1974 was among the first researchers to report these variations. He observed that fixation duration increases while saccade size decreases over the course of scene inspection. This early observation was confirmed by a number of studies [START_REF] Over | Coarse-to-fine eye movement strategy in visual search[END_REF][START_REF] Tatler | Systematic tendencies in scene viewing[END_REF]. The variation in the duration of visual fixations is contingent upon factors such as the quality of the stimulus and the number of objects in the scene (as explained in section 3.1). However, this variance might be explained by functional differences in the fixations. To investigate this point, Velichkovsky and colleagues [START_REF] Velichkovsky | Heterarchy of cognition: The depths and the highs of a framework for memory research[END_REF][START_REF] Unema | Time course of information processing during scene perception: The relationship between saccade amplitude and fixation duration[END_REF] conjointly analysed the fixation duration and the subsequent saccade amplitude. They found a non-linear distribution indicating that i) short fixations are associated with long saccades and, conversely, ii) longer fixations are associated with shorter saccades (figure 6 in [START_REF] Unema | Time course of information processing during scene perception: The relationship between saccade amplitude and fixation duration[END_REF]. This dichotomy permits us to disentangle focal-ambient fixations, using the terminology introduced by Trevarthen in 1968. Ambient processing is characterized by short fixations associated with long saccades. This mode might be used to extract contextual information in order to identify the whole scene. Focal processing is characterized by long fixations with short saccades. This mode may be related to recognition and conscious understanding processes. [START_REF] Pannasch | On the control of visual fixation durations in free viewing of complex images[END_REF] proposed the classification of fixations based on the amplitude of previous saccades. If the preceding saccade amplitude is greater than a threshold, the fixation is assumed to belong to the ambient visual processing mode. Otherwise, the fixation belongs to the focal mode. The authors chose a threshold equal to five degrees of visual angle. This choice is justified by the size of the parafoveal region in which visual acuity is good. Recently, [START_REF] Follet | New insights into ambient and focal visual fixations using an automatic classification algorithm[END_REF] proposed an automatic solution to classify visual fixations into focal and ambient groups. From this classification, they computed two saliency maps, one composed of focal fixations and the other based on ambient fixations. By comparing these maps to model-predicted saliency maps, they found that focal fixations are more bottom-up and more centred than ambient ones.

Conclusion

This paper provides an extensive overview of the different ways of analysing diachronic variables from eye-tracking data, because they are generally under-used by researchers. These diachronic indicators are scanpaths or saliency maps generated to represent the sequence of fixations over time. They are usually provided by eye-tracking software for illustrative purposes, but no real means to compare them are given. This paper aims to fill that gap by providing different methods of comparing diachronic variables and calculating relevant indices that might be used in experimental and applied environments. These diachronic variables give a more complete description of the visual attention time course than synchronic variables, and may inform us about the underlying cognitive processes. The ultimate step would be to relate the visual behaviour recorded with eye-trackers accurately to the concurrent thoughts of the user.

Despite looking at many analysis methods, some variables are still ignored (fixation duration, pupil diameter, etc) and it is very challenging to study the way these variables can be taken into account within diachronic data. A great improvement was recently made by combining eye movements with other techniques such as fMRI or EEG. For example, the development of EFRP (Eye-Fixation-Related Potentials) that tries to associate the displacement of the eye with some brain wave components [START_REF] Baccino | Eye movements and concurrent ERP's: EFRPs investigations in reading[END_REF][START_REF] Baccino | Eye-Fixation-Related potentials: Insight into parafoveal processing[END_REF]) is a first step in that direction. But other tracks should be explored, such as EDR (Electro Dermal Response) or ECG (electrocardiography). We are confident that researchers in this area will find new ways to go further in order to have a more complete understanding of human behaviour.

Figure 1 :

 1 Figure 1: Computation of a string edit distance to align the sequences ABCDE and ABAA recorded on a web page. First, AOIs are segmented and coded by letters (A, B, C…). Second, the substitution operations are carried out. The total cost is equal to 3 (the minimum number of substitutions), and normalized to the length of the longer string, here 5, yielding an edit distance between the two strings of d=(1-3/5) = 0.4.

Figure 3 :

 3 Figure 3: From left to right, (a) Original (b) Fixation map with red fixation points (c) Heat map (red spots represent the most visually salient areas of the picture) and (d) Saliency map.

  Figure 5 : thresholded saliency maps to keep the top percentage of salient areas. From left to right, 2, 5, 10 and 20%.

Figure 5

 5 Figure 5 gives an illustration of the thresholding operation on the Parrot picture (Figure 3). The first continuous saliency map (b) of Figure 4 is thresholded to keep 20 % of the image (T 20,G

Figure 7 .

 7 Figure 7. There are different methods to compute the AUC. The simplest ones are based on the left and right Riemann sums. The left Riemann sum is illustrated by Figure7. A more efficient approximation can be obtained by a trapezoid approximation: rather than computing the area of rectangles, the AUC is given by summing the area of trapezoids. In our example, the AUC value is 0.83.

Figure 6 :

 6 Figure 6 : Classification result (on the right) when considering a 20% thresholded ground truth (left picture) and a prediction (middle picture). Red areas are True Positives, green areas are False Negatives, and blue areas are False Positives. Other areas are True Negatives.

Figure 7 :

 7 Figure 7: Pseudo code to perform an ROC analysis between two maps (a), ROC curve (b) and the confusion matrix (c). The AUC is approximated here by a left Riemann sum as illustrated in (b).

Figure 8 :

 8 Figure 8: Example of ROC analysis. Red and green dots are the fixations of two observers for the Parrots image. These dots are drawn on a thresholded saliency map. On the left hand side the hit rate is 100% whereas the rate is 50% for the example on the right hand side.

Figure 9 :

 9 Figure 9: Distribution of the first five fixations for five observers, combined across 7 pictures (top of figure).

𝛼 2 ×

 2 100% of each tail of the distribution where α is the confidence level.

Figure 10 :

 10 Figure 10: The red dots correspond to the visual fixations of one observer viewing the Parrots image. Other colours correspond to fixations of the same observer but for three different pictures. Control fixations are chosen from this set of fixations.

Figure 11 :

 11 Figure 11: Example of NSS computation: the heat map is a normalized version of the model-predicted saliency map with a zero mean and unit standard deviation. A scanpath composed of 8 fixations (grey circles; the black one is the first fixation) is overlaid upon the standardized map. The normalized salience is extracted for each location. Values are shown in black next to the fixations.

Figure 13 :

 13 Figure 13: IOC as a function of the number of fixations for two different pictures. a) Error bars represent the standard error of the mean; b) examples of IOC extracted from pictures used by Le Meur, Baccino & Roumy (2011).

Figure 13

 13 Figure13(a) illustrates the IOC computation as a function of the numbers of fixations for two different pictures. The first picture (Parrots) represents two parrots that are visually salient. The second picture (Stream) is a mountain landscape where no element strongly attracts our attention. For both pictures, the congruency decreases over time as expected. The highest congruency, observed at the beginning of the viewing, is likely to be due to the bottom-up influences and the central bias. It is interesting to emphasize that the congruency for the Parrots image is significantly higher than that observed for the Stream image. This observation is mainly due to the attractiveness of the two contents. When there is nothing in the scene that catches our attention, observers are not 'unconsciously-constrained' and they just do not explore the visual scene in the same way. Figure13(b), extracted from Le Meur, Baccino & Roumy, 2011, shows other examples of IOC values.

  Figure13(a) illustrates the IOC computation as a function of the numbers of fixations for two different pictures. The first picture (Parrots) represents two parrots that are visually salient. The second picture (Stream) is a mountain landscape where no element strongly attracts our attention. For both pictures, the congruency decreases over time as expected. The highest congruency, observed at the beginning of the viewing, is likely to be due to the bottom-up influences and the central bias. It is interesting to emphasize that the congruency for the Parrots image is significantly higher than that observed for the Stream image. This observation is mainly due to the attractiveness of the two contents. When there is nothing in the scene that catches our attention, observers are not 'unconsciously-constrained' and they just do not explore the visual scene in the same way. Figure13(b), extracted from Le Meur, Baccino & Roumy, 2011, shows other examples of IOC values.

Figure 14 :

 14 Figure 14: Models performance tested on (a) the Le Meur dataset (top); (b) the Bruce dataset (bottom).All models perform better than chance and worse than humans. Judd's model gives the best performance on average.

7

  Limitation: Do visual fixations have the same meaning?

Table 1 :

 1 Eye-tracking datasets freely available on the web.

	Dataset name	Number of participants	Task	Number of pictures	Stimuli presentation time
	Le Meur	30	Free viewing	40	14s
	Kootstra	31	Free viewing	100	-
	Bruce	20	Free viewing	120	4s
	MIT	15	Free viewing	1003	3s
	DOVES	29	Free viewing + task	101 (greyscale)	5s
	FIFA	27	Free viewing + task	180	2s
	Ehinger	14	Task	912	-

Table 2 :

 2 NSS scores for four state-of-the-art saliency models on the Le Meur and Bruce datasets. SEM is the Standard Error of the Mean. A high average NSS value indicates a good prediction. Prediction of computational models may depend on the dataset. For instance, Itti's model performs much better on Bruce's dataset than on Le Meur's dataset. Judd's model gives similar results for both datasets.

	Model / Dataset	Le Meur	Bruce
	Itti	0.60+ /-0.10	0.99+ /-0.05
	Le Meur	0.77+ /-0.13	0.87+ /-0.03
	Bruce	0.60+ /-0.09	0.72+ /-0.04
	Judd	0.82+ /-0.11	0.87+ /-0.05


Alexandre Ganachaud
email: alexandre.ganachaud@ird.fr

Alexander Sen Gupta

Jaclyn N Brown

Karen Evans

Christophe Maes

Les C Muir

Felicity S Graham

A Sen Gupta

Projected changes in the tropical Pacific Ocean of importance to tuna fisheries

Keywords: Climate Change, Ocean Circulation, Stratification, Thermocline, Sea Surface Temperature, Warm Pool, Equatorial Upwelling, Tuna, Fisheries, Skipjack, CMIP3, CMIP5 PACS 92.70.Jw, 92.70.Gt, 92.70.Er, 92.10.Gk

Future physical and chemical changes to the ocean are likely to significantly affect the distribution and productivity of many marine species. Tuna are of particular importance in the tropical Pacific, as they contribute significantly to the livelihoods, food and economic security of island states. Changes in water properties and circulation will impact on tuna larval dispersal, preferred habitat distributions and the trophic systems that support tuna populations throughout the region. Using recent observations and ocean projections from the CMIP3 and preliminary results from CMIP5 climate models, we document the projected changes to ocean temperature, salinity, stratification and circulation most relevant to distributions of tuna. Under a business-as-usual emission scenario, projections indicate a surface intensified warming in the upper 400 m and a large expansion of the western Pacific Warm Pool,

Introduction

Tuna harvested from the tropical Pacific Ocean and the ecosystems that support them depend intimately on the oceanic environment. The distribution of these marine predators is linked to the horizontal displacement of water of a suitable temperature, and to vertical changes in the depth of the mixed layer. Areas where currents diverge or converge are also of major importance because they are associated with the thermal fronts, upwelling and eddies which enhance local productivity and, in turn create important foraging areas.

Tuna make up the overwhelming majority of fish harvested in the tropical Pacific Ocean and contribute significantly to the livelihoods, food and economic security of many island nations throughout the region (Bell et al. 2011, chapt. 1). A number of species are of particular importance including: skipjack (Katsuwonus pelamis, which comprise more than 60% of the catches in the region), yellowfin (Thunnus albacares), bigeye (Thunnus obesus) and albacore tuna (Thunnus alalunga). The response of these predators to oceanic conditions varies between species and life stage (i.e. larvae, juvenile and adult) and is primarily associated with the physiology of the species and how the species they prey on -the microneckton-are distributed. Of the four species, skipjack tuna appear to have the most limited physiology, with adults predominately occurring in waters with a temperature range of 20 • -29 • C and with a dissolved oxygen concentration above ∼ 4.0 ml/l. Temperature and oxygen tolerances are generally wider in the other three species with albacore predominantly occurring in cooler waters than yellowfin and bigeye, and bigeye able to tolerate lower dissolved oxygen concentrations than yellowfin and albacore [START_REF] Bushnell | Oxygen transport and cardiovascular responses in skipjack tuna (Katsuwonus pelamis) and yellowfin tuna (Thunnus albacares) exposed to acute hypoxia[END_REF]. Tuna are most sensitive to water temperatures during their larval stages, with thermal tolerances widening as they grow. The distribution of primary production is also important for earlier life stages of tuna and micronekton. This distribution will depend on both the physical and chemical properties of the ocean and in particular, the availability of oceanic nutrients. These are abundant at depths, most often below the biologically-active layer (euphotic zone), and the supply to this euphotic zone depends upon oceanic currents and mixing processes.

Five ecological areas will guide our analysis, each with distinct vertical hydrological structure and ecosystems (Figure 1a): the Pacific Equatorial divergence (PEQD), the Pacific Warm Pool (WP), the North Pacific Tropical Gyre (NPTG), the South Pacific Subtropical Gyre (SPSG) and the Archipelagic Deep Basins (ARCH; [START_REF] Longhurst | Ecological geography of the Sea[END_REF]. The time-dependent boundaries between each of these provinces are predominantly associated with the convergence zones of ocean surface currents. Of particular relevance to tuna catches in the western Pacific is the boundary that exists between WP and PEQD, the position of which is influenced by the El Niño Southern Oscillation (ENSO), moving eastward during El Niño and westward during La Niña by thousands of kilometres [START_REF] Maes | Characteristics of the convergence zone at the eastern edge of the pacific warm pool[END_REF]). More than 90% of skipjack tuna catches are located within the WP and demonstrate associated shifts in distributions with ENSO (Lehodey et al. 1997, this issue). Any future changes to the WP or to the characteristics of ENSO could therefore have a large influence on the availability of tuna stocks to Pacific island states and to their associated economies.

This paper provides fisheries scientists and managers with information on future projections of the physical properties of the tropical Pacific. We focus on the parameters that are of relevance to tuna, and in particular, skipjack tuna, as this species comprises the greatest proportion of tuna catches throughout this region [START_REF] Lehodey | El Niño Southern Oscillation and tuna in the western pacific[END_REF]. We first introduce the climate models used to project the ocean state (Section 2), then examine changes to surface and subsurface temperatures, vertical stratification, major ocean currents and upwelling, as well as issues associated with future projection of ENSO (Section 3). Finally, we describe some of the limitations inherent in those projections and discuss the implications for tuna as a result of the projected physical changes (Section 4).

Methods

To examine future changes to the ocean we use output from the third Coupled Model Intercomparison Project (CMIP3), archived by the Program for Climate Model Diagnoses Intercomparison (PCMDI). The CMIP3 repository contains 24 nominally independent models encompassing a broad range of resolutions and a variety of physical parameterizations (Sen Gupta et al. 2009, Table 1). Ocean horizontal resolutions vary from 1

4

• to 5 • . Unless otherwise stated we compare present day conditions from the end of the 20th century with future projections for 2100 and from the business as usual SRESA2 emissions scenario [START_REF] Nakicenovic | IPCC special report on emissions scenarios[END_REF]. Projections for 2100 are calculated by taking a linear trend in a variable over the full 21st century-unless otherwise stated, to minimize the effect of low frequency variability that can distort the climate change signal. The significance of projected changes is established using a standard t-test. This assumes independence among the models, which may lead to overestimation of the areas where the changes should be considered significant (Pennel and Reichler 2011). Some climate models exhibit spurious trends, unrelated to external forcing, that may affect the estimated climate change signal. This drift is of a different sign, magnitude and spatial distribution in different models and as such has a negligible effect on averages across multiple models (Sen Gupta et al 2012b;[START_REF] Tebaldi | The use of the multi-model ensemble in probabilistic climate projections[END_REF]. As the analysis below generally focusses on the multi-model mean (MMM), we have neglected errors introduced by model drift. However, when considering individual models, drift is likely to be a small additional source of uncertainty.

We have adopted the approach of examining the largest number of models possible for each diagnostic. As such, the number of models used may vary in a given analysis (Online Table S1). Given the importance of ENSO behaviour on tuna, we also include preliminary analysis on the currently available subset of models taking part in CMIP5 when discussing projected changes to interannual variability.

Results

Sea surface temperatures

Under typical conditions, trade winds push the warmest tropical waters to the western side of the Pacific basin, forming the WP, where sea surface temperature (SST) generally exceeds 29 • C (Fig. 1a). Along the equator and the coast of South America, the prevailing winds cause equatorial and coastal upwelling, respectively, forming the PEQD, a relatively cool tongue of water (SST< 20 • C) that extends across the Pacific and meets the eastern edge of the WP. At 10 • N, a band of warm water extends across the entire basin, below the intertropical convergence zone (ITCZ). At 10 • S, a similar band only extends to 140 • W a result of the South Pacific Convergence Zone (SPCZ) being confined to the western Pacific. SST varies seasonally by ±2 -3 • C near 25 • N/S (e.g., Bell et al. 2011, chapt. 2), but by only ±1 • C within 15 • of the equator, where instead interannual variability dominates (±2 • C, Section 3.5). Over the 1955-2003 period, warming trends of about 0.2 • to 1 • C have been observed in the western Pacific [START_REF] Cravatte | Observed freshening and warming of the western pacific warm pool[END_REF], thereby increasing the horizontal extent of the WP waters with temperatures greater than 29 • C.

The CMIP3 MMM reproduces the observed largescale temperature field (Fig. 1b). However, there is a systematic cold offset, or bias in the central and western equatorial Pacific (i.e. the bias is present in almost all the models) with simulated temperatures up to 2 • C too low in the MMM, with more extreme biases in individual models (Fig. 1c). Conversely, on the eastern boundaries, upwelling is too weak and temperatures too warm.

The projections for 2100 suggest a statistically significant MMM warming that exceeds 2.5 • C over large parts of the basin. The physical mechanisms driving the spatial pattern of this warming are reasonably well understood: a projected increase in the south easterly trade winds (south of the equator), leads to a relatively slow warming rate in the southeastern Pacific [START_REF] Xie | Global warming pattern formation: Sea surface temperature and rainfall[END_REF] while the enhanced equatorial warming is related in part to a weaker Walker circulation [START_REF] Dinezio | Climate response of the equatorial pacific to global warming[END_REF].

The WP region is projected to warm in all models, with the warmest 9 million square km (which correspond to the area of water enclosed by the 29 • C isotherm in observations) increasing by ∼ 2.7 ± 0.4 • C (±0.4 • is the standard deviation among the CMIP3 models). This is associated with the eastern edge of the WP (as defined by the 29 • C isotherm) moving east by 6000 ± 2000 km (Fig. 2a). Projections also show that the warming will be greatest in the cold tongue (> 3 • C east of 140 • E, Section 4).

It is often more pertinent to define the WP edge by a certain salinity surface (isohaline) than an isotherm (i.e. 29 • C). This alternative definition better reflects the position of the front generated by strong rainfall in the west and is a better proxy for the convergence zone between WP (or fresh pool) and PEQD. The WP surface salinity has decreased substantially over the past 50 years (1955-2003, 0.34 PSU, Cravatte et al. 2009), although the extent to which natural variability has contributed to this is uncertain. While the CMIP3 models also generally simulate a freshening in the WP over this period, the magnitudes of change are considerably smaller. The CMIP3 models generally project further freshening (Online Figure S1) due to increased rainfall over the WP and convergence zones. Projections suggest that the edge of the fresh pool will be displaced eastward by ∼ 500 -2000 km by 2100, much less than the WP (Fig. 2a).

Vertical temperature and salinity structure

Near the equator, WP waters extend to depth of ∼ 100 m, overlying a steep thermocline situated between 100 and 250 m (Fig. 3a). In PEQD the mixed layer and the thermocline are considerably shallower. Observational data from 1950 to 2008 show that surface warming of the central Pacific has penetrated vertically from 50 m to 100 m at the equator to the top of the thermocline (> 1 • C down to 100 m, Fig. 3c). Surprisingly, in the western Pacific, water temperatures have actually decreased by up to 2 • C through the thermocline depths (100 -150 m, Fig. 3c). It has been suggested that this cooling is related to a weakening of the easterly equatorial winds that through dynamical ocean processes cause an adiabatic lifting of the thermocline [START_REF] Han | Interpretation of tropical thermocline cooling in the indian and pacific oceans during recent decades[END_REF]. [START_REF] Vecchi | Examining the tropical pacific's response to global warming[END_REF] note however that the brevity of the data record precludes attributing the changes in the Walker circulation to anthropogenic Global Warming. The depth of the observed warming varies with latitude (Online Figure S2a), reaching 400 to 600 m in the northern hemisphere subtropics (not shown; see Bell et al. 2011, chapt. 3). The subsurface cooling noted in the western equatorial region (10 • S-10 • N) is evident between ∼ 100 -250 m. It extends southwards reaching depth exceeding 800 m in the southern hemisphere subtropical gyre (not shown). This enhanced vertical contrast is also observed in the salinity change, with a decrease in the WP and an increase in the thermocline waters below the WP [START_REF] Durack | Fifty-year trends in global ocean salinities and their relationship to broadscale warming[END_REF].

In the equatorial band, projected MMM warming is generally confined within the model mixed layer, reaching +2.5 • C at 50 m and ∼ +1 • C at 100 m-which is of a similar magnitude to observed trends of the 20th century (Fig. 3c,4a). Poleward of 10 • N/S, the projected warming penetrates to considerably greater depths than at the equator, with a warming of ∼ 1.5 • C at 200 m (Fig. 2e,4c). Below 100 m, and west of about 140 • W, very weak warming or moderate cooling is projected to occur straddling the equator (Fig. 3e), at about 160 m and within 10 • N/S. The cooling pattern is also similar to the observed trends and is again consistent with a shoaling of the thermocline driven by a weakening of the equatorial trade winds. Projected salinities also show an enhanced freshening in the WP (Online Figure S1) down to ∼50 m (not shown).

Stratification and mixed layer depth

Stratification, which is associated with the density difference between the surface and the deep oceanic layers, characterizes the stability of the water column, and therefore the potential extent of vertical exchange of properties such as nutrients or oxygen. Its interplay with the mixed layer, the region of homogenized temperature, salinity, oxygen and nutrients that extends from the surface to the mixed layer depth (MLD), has important influences on primary production (Section 4).

The combined effect of the surface intensified warming and, in the WP, decrease in salinity, has lead to a reduction in water density in the upper layers, consequently increasing ocean stratification over the upper 200 m [START_REF] Cravatte | Observed freshening and warming of the western pacific warm pool[END_REF]). In the MMM the stratification tends to be too strong in the southern hemi-sphere, possibly attributable to the tendency for the SPCZ to extend too far east (Fig. 5a-c). Under the ITCZ in the north, the MMM does not capture the full strength of the stratification nor is it well represented in the regions of the cold tongue bias.

Stratification is projected to continue to increase across most of the tropical Pacific (Fig. 5d), with the largest increases occurring in a wide area around the WP and under the SPCZ and ITCZ where both temperature and salinity effects reinforce each other. This increase in stratification is projected to occur at depths of between 20 and 150 m in the WP and at 20-100 m in PEQD depending on location (thin lines on Fig. 4; see also Online Figure S3). In the SPSG and ARCH, the depth range of increase widen to between 30 and 450 m, while increases in the NPTG occur predominantly between the surface and 250 m (Fig. 4c; Online Figure S4).

The observed annual maximum depth of the mixed layer (MMLD) does not generally exceed 100 -120 m in the highly stratified WP (Online Figure S5a), but reaches considerably deeper in the centre of the NPTG and SPSG, where high-salinity waters are cooled during the winter, causing convection and much deeper mixing. The MMLD is projected to shoal in most of the tropical Pacific (Online Figure S5d), with the most pronounced shoaling occurring in the WP (∼ 10 -25 m), and in the western half of both tropical gyre provinces poleward of ∼ 15 • N/S (> 40 m). Under the SPCZ and ITCZ, where model bias is highest (Online Figure S5c), projections are not consistent. Projected increases in the winds near 25 • S / 120 • W will result in the MMLD increasing by up to 40 m in this region.

Ocean currents and upwelling

Currents within the upper water column are expected to change across the tropical Pacific in the future, particularly as a result of weakened equatorial and northeasterly trade winds and strengthened south easterly trade winds [START_REF] Bell | Vulnerability of Tropical Pacific Fisheries and Aquaculture to Climate Change, secretariat of the pacific community edn[END_REF]; Sen Gupta et al 2012b, chapt 3). The MMM projections for large-scale surface currents (0-50 m, Fig. 6b) show relatively small absolute changes away from the equatorial band (12 • S-10 • N). Near 9 • S the upper South Equatorial Counter-Current (SECC), which flows eastward across the Pacific, is projected to decrease substantially. Its surface flow is projected to turn towards the south, thereby reducing its penetration to the east. The South Equatorial Current (SEC), which flows westward across the Pacific, is projected to decrease in strength dramatically within 3 • of the equator, from ∼ 30 cm.s -1 to ∼ 20 cm.s -1 in 2100. Just underneath, the upper part of the eastward flowing Equatorial Under-Current (EUC) is projected to increase substantially, leading to a ∼ 20 m shoaling of the EUC core (Fig. 6c). At ∼ 7 • N, a decrease in the upper part of the eastward flowing North Equatorial Counter-Current (NECC) is projected, while the broad North Equatorial Current (NEC) demonstrates no significant changes. The location of the fronts that lie between the westward equatorial SEC, the eastward SECC and NECC, which define the poleward boundaries of PEQD (Fig. 6b), also demonstrates little projected change.

Upwelling at the equator is projected to decrease in intensity, as are the neighboring two regions of downwelling near 4 • N/S (Figure 7). This results in little projected change in the net upwelling integrated between 9 • S and 9 • N, which is most important to equatorial temperature variations [START_REF] Zhang | Decadal variability of the shallow pacific meridional overturning circulation: Relation to tropical sea surface temperatures in observations and climate change models[END_REF].

Interannual variability

We have explicitly attempted to isolate the global warming signal and remove the effects of year to year variability in our analyses. However, long-term trends are likely to impact ENSO characteristics. At present, there is little consistency across the models with regard to possible changes in the magnitude or frequency of ENSO variability. This is in part related to the significant biases in the climate models and their representation of ENSO (Brown et al this issue). It also relates to the fact that the characteristics of ENSO depend upon a number of competing feedbacks that act to alter its amplitude [START_REF] Collins | The impact of global warming on the tropical pacific ocean and El Niño[END_REF]. Both CMIP3 and preliminary results from CMIP5 models demonstrate both increases and decreases in frequency and amplitude of ENSO, generally within 20% of their 20th century values (Fig. 2c).

Observations show that the depth of the thermocline decreases during El Niño events and increases during La Niña events in the west Pacific, with the opposite occurring in the east Pacific, on the order of 20 -30 m (Online Figure S6). Average shoaling of the thermocline is projected to be in the order of ∼ 10 -20 m, by 2100. This means that during El Nino events the thermocline could be ∼ 40 m shallower than the 20th century mean position in the WP.

As noted above, the mean longitude of the eastern edge of the warm/fresh pool, which is subject to large interannual excursions due to ENSO, will move eastwards in the future. While the projected strength of ENSO (Figure 2c) either increases or decreases depending on the model, projected variability in the location of the fresh pool edge shows either little change or an increased variability (Figure 2b). This highlights non-linearities that are still to be understood in model projections of ENSO.

Discussion

Significance and limitations of this analysis

There are a number of known issues when examining the future projections of the climate of the tropical Pacific. Primary amongst these is the cold tongue bias: most of the climate models are too cold along the equator and the cold tongue extends too far to the west. This error also manifests in the simulation of precipitation and ENSO, with impacts on the WP structure. Another major bias is associated with the SPCZ; it is zonally elongated and extends too far to the east in most models [START_REF] Brown | Evaluation of the south pacific convergence zone in IPCC AR4 climate model simulations of the twentieth century[END_REF]. As a result, precipitation and wind patterns are distorted, with consequences for ocean salinity and currents in the southwest Pacific. These biases induce offsets that must be considered when assessing future projections (Brown et al this issue).

Here, we have used output from the full suite of climate models taking part in CMIP3, to examine projections in key physical parameters that are known to affect the distribution of tuna species and their prey. We have also discussed some of the systematic biases that are evident in these properties and examined both historical and projected changes in tropical Pacific over the 21st century. The projections that we report should therefore be interpreted in relation to their location with respect to ecological provinces and dynamical features. For example, as the MMM cold tongue extends too far to the west, we expect that the large warming projected between 140 • E and of 160 • W may be exaggerated, and a more realistic projection might confine the largest warming to the east of 160 • W (Brown et al this issue).

While part of the observational trends discussed herein may be due to natural variability [START_REF] Vecchi | Examining the tropical pacific's response to global warming[END_REF], overall agreement with respect to the patterns of change with climate models forced with increasing greenhouse gas concentrations imply that the observed trends are in general consistent with scenarios of anthropogenic Global Warming. By averaging over multiple climate models, we are able to reduce non-systematic model biases and identify regions of model consensus. Averaging, by its nature, removes the effect of outliers. More (or less) extreme projected climate changes are nonetheless feasible. As an example, while the projected MMM SST warming in the warm pool area (defined as the warmest 9.10 6 km 2 ) is 2.7 • C, three models show warming that is over 3.2 • C.

Relevance of projections for tuna fisheries

The extension and warming of the WP and associated stratification of the water column is likely to influence tuna distributions and productivity across the tropical Pacific region both directly, through changes in preferred thermal and dissolved oxygen habitats and indirectly, via changing distributions and abundances of tuna forage. Overall, it is anticipated that tuna populations will re-distribute across the tropics and subtropics and overall productivity throughout the region will decrease (Lehodey et al. this issue). However, climate models can show significant inter-model differences in the present day fidelity and their projections of future mean climate and variability and contain a number of significant biases.

Skipjack tuna have a preferred thermal range of 20 to 29 • C, with few tuna observed at temperatures beyond 17-30 • C. Along the equator, the boundary of upper preferred temperatures (T= 29 • C) is projected to move eastward by 6000 ± 2000 km (Fig. 8a), accounting for the large inter-model spread. Below the surface, this isotherm is projected to move down to ∼ 100 m west of 170 • W, while the depth of the lower preferred temperature (T= 20 • C) shows little change. In the eastern part of the basin, the 20 • C isotherm is projected to deepen slightly, extending waters of a suitable thermal range to skipjack. In the west Pacific, surface waters warmer than 29 • C are projected to expand to 20 • S and more than 20 • N (Fig. 8b, see also Online Figure S7). Assuming similar thermal tolerances into the future and on the basis of preferred thermal habitats alone, skipjack tuna are likely to demonstrate substantial eastward and poleward displacement in their habitat, consistent with existing tuna modeling projections (Lehodey et al. this issue).

The NPTG, SPSG and ARCH provinces are characterized by a thick thermocline and a MMLD which is situated below the euphotic zone Consequently these regions are characterised by low levels of productivity. Projected increases in stratification of the water column in these provinces will increase dynamical barriers, further limiting exchange between the euphotic layer, and the deep, nutrient-rich layer (the nutricline), thereby limiting further primary productivity. Conversely, in the WP, the MMLD has little seasonality and is situated within the euphotic zone so that primary productivity can take place beneath it. Projected shoaling of the MMLD may therefore increase nutrient exposure to the sun, leading to enhanced primary production. These processes however, will ultimately depend on the future distribution of nutrients and the position of the nutricline, whose estimation necessitates full biogeochemical models. Recent efforts coupling a biogeochemical model with one of the CMIP3 models projected a deepening in the nutricline, below the WP, resulting in overall lowering of primary production (Bell et al. 2011, chapt. 4).

Most of the models are able to simulate the major current systems in the tropical Pacific. There is however large inter-model differences in the strength and location of the currents. The SECC in particular is poorly represented in many of the models, and non-existent in some of them as a consequence of the SPCZ biases, discussed above. Projected changes in the strength of the major currents across the Pacific will affect biological activities and associated fisheries by influencing horizontal transport of nutrients essential to primary productivity, mixing processes supplying nutrients to the euphotic zone, as well as influencing larval dispersion and connectivity [START_REF] Munday | Climate change and coral reef connectivity[END_REF]. The largest projected changes to the upper part (0 -50 m) of the major currents are those close to the equator where the SECC, NECC and the equatorial part of the SEC are projected to weaken. Sub-surface the upper EUC is projected to shoal substantially and intensify in the western and central Pacific. As the EUC is the primary conduit of iron to the equatorial region [START_REF] Mackey | Iron in the western pacific: a riverine or hydrothermal source for iron in the equatorial undercurrent?[END_REF], such change could increase the productivity in the eastern and central Pacific surface layers, where iron appears to be the limiting nutrient (Sen Gupta et al 2012a).

Dissolved oxygen concentrations are also important in determining future projections of tuna distributions, but are not available in the CMIP3 models. This thereby limits the ability to examine interactions between temperature and oxygen concentrations which may influence future tuna habitats.

Conclusions

A number of robust physical changes are projected for the tropical Pacific which are likely to influence tuna distribution and associated fisheries. While there are significant biases in surface temperature in the models which may affect the future distribution of temperature change, substantial warming is evident everywhere, except in the south eastern basin. This places much of the western tropical Pacific outside of the temperature range at which skipjack tuna exist today, reducing dramatically their preferred habitat at low latitudes, and also increasing waters within their thermal range at higher latitudes.

While we find some robust changes in the meanstate, tuna distributions are also subject to dramatic interannual variability and therefore likely to be very sen-sitive to any future changes in ENSO. Significant biases exist in the simulation of ENSO in the CMIP3 models and there is little consistency with regards to what these models suggest might happen to ENSO characteristics in the future. Considerable work has gone into improving the simulation of ENSO in the new CMIP5 models. Nevertheless, based on our preliminary results from a subset of these models we are still no closer to obtaining a model consensus.

Future warming will be more intense in the surface ocean layers and, along with significant freshening driven by increased rainfall in the west, will lead to an increase in stratification. Increased stratification and weaker winds in the equatorial region is projected to cause a shoaling of the mixed layed in most areas except the southeastern Pacific. This will likely impact both primary productivity and associated tuna forage: In the subtropics, nutrient availability is likely to be reduced. Under the WP, uncertainty in the future depth of the nutricline and associated specific processes make projections of productivity unclear. In the PEQD, the EUC shoaling has a potential to increase productivity; however the combined effect of increases in stratification and reduction in upwelling are not clear.

A number of models participating in CMIP5, will include biogeochemical processes. This will allow direct projections for both physical parameters and nutrient concentrations, oxygen and phytoplankton. Assessing the fidelity of these biogeochemical parameters will be a challenge, however, given the relative sparsity of chemical and biological data compared to physical data such as temperature. This highlights the urgency of initiating and maintaining new continuous monitoring programs, particularly in the South Pacific where there are no long term biogeochemical monitoring programs Increased Frequency

Decreased Frequency

Increased Amplitude Decreased Amplitude . Amplitude corresponds to one standard deviation. Frequency was computed as a weighted average of NINO3.4 frequency by spectral power for all frequencies whose power was significant at 95% level. Panels (b) and (c) contain addition preliminary results from the CMIP5 archive. CMIP5 data is still in progress and modelling groups may replace certain data in the future. The shaded area delimits the 20% change.

CSIRO-Mk3-6-0 CanESM2 GISS-E2-R MIROC5 MRI-CGCM3 inmcm4 CSIRO-Mk3-6-0 CanESM2 GISS-E2-R MIROC5 MRI-CGCM3 inmcm4

Fig. 1

 1 Fig. 1 (a) Observed climatological SST from CARS with warm pool boundary (contour at 29 • C) and oceanic provinces superimposed, (b) multi-model mean (MMM; 1980-1999) SST, based on 19 models, (c) MMM bias compared to observations, with a strong cold bias along the equator (d) change in SST per 100 years for 21st century of the simulation under the SRES A2 scenario (based on a linear trend line). Quantitative model agreement is ubiquitous at 95% based on a Student t-test; Mottling indicates regions where there is at least 80% agreement on the sign of the temperature change when the area mean change has been removed (i.e. a measure of the agreement in the spatial pattern of change). Oceanic provinces are as follows: Pacific Equatorial divergence (PEQD), Pacific Warm Pool (WP), North Pacific Tropical Gyre (NPTG), South Pacific Subtropical Gyre (SPSG) and the Archipelagic Deep Basins (ARCH).

  Fig.2Changes in WP and ENSO variability. a) change in the longitude of the equatorial (3 • S to 3 • N) warm pool and fresh pool edge over the 21st century in individual CMIP3 models. b) change in the variability of the fresh pool edge ( • longitude, 1 standard deviation). The fresh pool and warm pool are defined by model specific isohalines and isotherms, respectively, chosen such that the area of bounded freshwater or warm water was equivalent to the area contained within the T= 29 • C isotherm and S= 34.8 isohalines in the observations from CARS. c) Ratio of 2050-2100 to 1950-2000 frequency and amplitude of temperature anomalies (NINO3.4 index). Amplitude corresponds to one standard deviation. Frequency was computed as a weighted average of NINO3.4 frequency by spectral power for all frequencies whose power was significant at 95% level. Panels (b) and (c) contain addition preliminary results from the CMIP5 archive. CMIP5 data is still in progress and modelling groups may replace certain data in the future. The shaded area delimits the 20% change.

Fig. 3

 3 Fig. 3 Equatorial potential temperature averaged 5 • S to 5 • N ( • C) as a function of longitude and depth. a) Observed longterm mean b) Model bias (MMM-observations). c) Observed linear temperature trend over 50 years. d) MMM trend over 1950-2000. e) Projected trend in temperatures 2000-2100 (color shading). Overlaid (line contours) is the MMM 1950-2000 mean (a,b,d,e) and observed long-term mean (c). Observations are based on a new dataset (Durack and Wijffels 2010; Johnson and Wijffels 2011) that has used observations from 1950-2008 to compute the mean state and trends (expressed as a 1950-2000 trend). Multi-parametric regression technique was employed to reduce aliasing of the trend by various modes of climate variability including ENSO and the PDO.

Fig. 4

 4 Fig. 4 Average vertical temperature profiles for typical sites of each province ( • C). The thick grey line corresponds to the modern climatological temperature from CARS; the grey dashed line corresponds to projected MMM 100-year trend added to CARS (1 standard deviation model spread: shaded region). The MMM stratification, thin black lines, is indicated by Brunt-Väisälä frequencies, N 2 (i.e the measure of the strength of the vertical density gradient). Projected N 2 , dashed black line, is thickened where the change is significant at 95% level. Provinces are defined as follows a) WP 5 • S-5 • N / 155 • E-160 • E; b) PEQD 3 • S-3 • S / 140 • W-130 • W; c) SPSG 25 • S-20 • S 160 • W-140 • W.

Fig. 5

 5 Fig. 5 Same as Fig. 1 but for density difference between 200m and surface (kg.m -3 ) (a) observed (WOA09), (b) MMM, (c) MMM bias, (d) MMM projected change over the 21st century for A2 scenario. Mottling indicates where there is at least 80% agreement on the sign of the change when the area mean change has been removed (i.e. a measure of the agreement in the spatial pattern of change).

Fig

  Fig.7a) Vertical velocities of water (10 -6 m.s -1 ) at a depth of 50 m, for 1990(mean 1980 -1999), averaged between 170 • E and 110 • W, where a positive velocity corresponds to upwelling and a negative velocity to downwelling; and b) Climate driven changes to vertical velocities for 2100 (full line if significant at 95%; dashed line otherwise). Values for individual CMIP3 models are indicated in the maxima (small circles). The strong equatorial upwelling shows decreases that are immediately balanced by strong downwelling near 4 • N-4 • S so that net upwelling between 9 • S and 9 • N does not change.

Acknowledgements Most of the analyses were based on the international modeling groups participating in IPCC AR4 and AR5, the Program for Climate Model Diagnosis and intercomparison (PCMDI). The IPCC Data Archive: allowed efficient analysis of the massive numerical model outputs. Much of the work described in this paper was commissioned by the Secretariat of the Pacific Community to help assess the vulnerability of fisheries and aquaculture in the tropical Pacific to climate change. J. Lefèvre's technical help was most valuable. We thank Johann Bell for his guidance on how to present the analyses to make them easy for fisheries scientists to use. K. Richard and N. Holbrook also provided valuable feedback at early stages of this work.

This work was partly funded by ANR project ANR-09-BLAN-0233-01; it is a contribution to the CLIVAR/SPICE. Part of the research was conducted with the support of the Pacific Australia Climate Change Science and Adaptation Program (PACCSAP), a program managed by the Department of Climate Change and Energy Efficiency in collaboration with AusAID, and delivered by the Bureau of Meteorology and the Commonwealth Scientific and Industrial Research Organisation (CSIRO).

A. Ganachaud (corresponding author),