Functional sliced inverse regression to infer temperature , water vapour and ozone from IASI data

Technical note: Functional sliced inverse regression to infer temperature, water vapour and ozone from IASI data U. Amato, A. Antoniadis, I. De Feis, G. Masiello, M. Matricardi, and C. Serio Istituto per le Applicazioni del Calcolo “Mauro Picone” CNR, Napoli, Italy Laboratoire Jean Kuntzmann, Université Joseph Fourier, Grenoble, France Dipartimento di Ingegneria e Fisica dell’Ambiente, Università della Basilicata, Potenza, Italy European Centre for Medium-Range Weather Forecasts (ECMWF), Reading, UK Received: 7 January 2009 – Accepted: 10 March 2009 – Published: 23 March 2009 Correspondence to: I. De Feis (i.defeis@iac.cnr.it) Published by Copernicus Publications on behalf of the European Geosciences Union.


Introduction
The development of satellite high-spectral resolution infrared spectrometers is expected to improve quality and density of retrieval of atmospheric parameters.Both Numerical Weather Prediction and Earth's monitoring are expected to benefit from these new modern sensors.
The pioneer of this new generation of instruments has been the Interferometric Monitor for Greenhouse (IMG) gases sensor (Kobayashi et al., 1999), that flew on the ADvanced Earth Observing Satellite (ADEOS) platform from August 1996 to June 1997.At present, high-resolution infrared sensors on operational meteorological polar orbiters include the Atmospheric Infrared Sounder (AIRS) on the second Earth Observing System (EOS) polar orbiting platform, EOS-Aqua, launched in April 2002 (Aumann and Pagano, Correspondence to: I. De Feis (i.defeis@iac.cnr.it)1994), and the Tropospheric Emission Spectrometer (TES) on the AURA satellite launched in 2004 (Beer et al., 2001).In March 2002, the European Space Agency launched Envisat, an advanced polar-orbiting Earth observation satellite which provides measurements of the atmosphere, ocean, land, and ice.It is equipped with the Michelson Interferometer for Passive Atmospheric Sounding (MIPAS), that is a Fourier transform spectrometer for the measurement of highresolution gaseous emission spectra at the Earth's limb and operates in the near to mid infrared (Fischer et al, 1996(Fischer et al, , 2000)).
The new arrived in the family of high resolution infrared sensors is the Infrared Atmospheric Sounding Interferometer (IASI) (EUMETSAT, 1998) on the first European Meteorological Operational Satellite (METOP/1) launched in 19 October 2006.IASI is a Fourier Transform Spectrometer based on a Michelson Interferometer coupled to an integrated imaging system that observes and measures infrared radiation emitted from the Earth in the spectral range 3.62-15.5µm (645-2760 cm −1 ), covering the peak of the thermal infrared and particularly the intense CO 2 band around 666 cm −1 , with an apodized resolution of 0.5 cm −1 and a spectral sampling of 0.25 cm −1 .IASI characteristics have been specified to get observations, which are compatible in terms of sampling, resolution, accuracy and overall performances with the mission objectives of providing improved information on temperature, water vapour, ozone, cloud top pressure and temperature, cloud cover and cloud optical properties.
The large amount of observed data needs algorithms for the radiative transfer equation and its inversion specifically designed for IASI.In this paper we describe a statistical regression methodology for temperature, water vapour and ozone, (T , q, o), which exploits IASI observations.The Published by Copernicus Publications on behalf of the European Geosciences Union.
U. Amato et al.: Dimension-reduction for IASI radiances methodology is based on a suitable statistical dimension reduction technique, FSIR (Functional Sliced Inverse Regressions.FSIR (Amato et al., 2006) generalizes the well known Principal Component Analysis (PCA) (Jolliffe, 2002, e.g.) or Empirical Orthogonal Function (EOF) approach (see the recent review on EOF regression methods by Serio et al., 2009) and allows one to deal with functional models.Functional data analysis is about the analysis of information on curves or functions, and the radiance spectrum is a curve evaluated at fixed points that depend on the design of the interferometer.
FSIR needs to be trained on a suitable set of pairs: radiances, profiles.We have selected the profiles from the well known ECMWF (European Centre for Mediumrange Weather Forecasts) Chevallier data base (Chevallier, 2001).
The Chevallier data base documentation can be downloaded from the EUropean Organisation for the Exploitation of METeorological SATellites (EUMET-SAT) Satellite Application Facilities (SAF) Scientific Report web page: http://www.eumetsat.int/groups/pps/documents/document/002197.pdf;also proper documentation and the data-base as well can be obtained directly by F. Chevallier at f.chevallier@ecmwf.int.IASI synthetic radiances have been computed using the radiative transfer code σ -IASI (Amato et al., 2002).
The code σ -IASI has been extensively validated with the use of the NAST-I instrument (Cousins and Gazarick, 1999, e.g.).NAST-I is the National Polar-orbiting Operational Environmental Satellite System (NPOESS) Airborne Sounder Testbed Fourier Transform spectrometer, flying onboard the NASA aircrafts, ER-2 and Proteus.An extensive retrieval exercise with σ -IASI has been performed for the CAMEX/3 experiment (Convection and Moisture Experiment 3) (Carissimo et al., 2005(Carissimo et al., , 2006, e.g.), e.g.).More recently, σ -IASI has also been used within the EAQUATE campaign (European AQUA Thermodynamic Experiment), for analysis of NAST-I spectra (Proteus Flight) recoded over the Mediterranean Sea (Taylor et al., 2008;Grieco et al., 2007).The forward module σ -IASI has been also validated with the use of AIRS (Atmospheric Infrared Radiometer Sounder) data, flying onboard the Aqua satellite (Saunders et al., 2007).
The application of FSIR to IASI data has been exemplified through a series of IASI soundings recorded over the tropical basin.These observations have been FSIR-regressed for (T , q, o).To simplify the comparison of retrieval products with truth data, only clear-sky, sea-surface IASI soundings have been analyzed in this work.Truth data have been derived from the ECMWF analysis for the same date and location as the IASI soundings.The difference (retrieval-ECMWF) for temperature, water vapour and ozone has been evaluated, which has allowed us to assess the retrieval performance of FSIR.The paper also provides a comparison with a conventional EOF or PCA regression scheme (see also Serio et al., 2009).
The paper is organized as follows.Section 2 will deal with the mathematical aspects of FSIR.Application of FSIR to IASI data will be described and discussed in Sect. 4. Finally conclusions will be drawn in Sect. 5.

The regression model
Functional Sliced Inverse Regression (FSIR) (Amato et al., 2006) is a statistical tool to reduce the dimensionality based on the Sliced Inverse Regression (SIR) (Li, 1991), that permits to deal with functional models.It generalizes the Principal Component Analysis (PCA) using inverse regression.
In functional regression problems, one predicts a response variable Y from a set of variables R 1 , . . ., R d that are discretizations of a same curve R at points ν 1 , . . ., ν d , that is R j = R(ν j ), j = 1, . . ., d, where the discretization points ν j lie in some interval.In our case Y is the geophysical variable to be retrieved (surface temperature, temperature or gas concentration in a fixed layer) and R 1 , . . ., R d are the measurements, i.e., the radiances at wavenumbers ν 1 , . . ., ν d .Let us consider the model where m is a smooth link function of the K-dimensional Euclidean space, E K , into the one-dimensional space or real axis E 1 , that indicates what kind of relation exists between the geophysical variable to retrieve and the radiance and is assumed to be linear; ε is noise assumed independent of the curve R(ν); {β i (ν), i=1, . . ., K} are K orthonormal functions with respect to the usual inner product defined in the space L 2 ([ν min , ν max ]) and denoted with the symbol < •, • >.
Let R be the covariance operator of R(ν) and e the covariance operator of the conditional expected value, E(R|Y ), of the radiance given the geophysical variable to retrieve.
Provided that K<d, with d denoting the number of points in which the radiance curve is known, Eq. (1) says that the regression function depends on R(ν) only through K linear functionals of the explanatory process R(ν).Hence, to explain the dependent variable Y , the space of d explanatory variables can be reduced to a space with a smaller dimension K.The dimension reduction methods aim at finding the dimension K of the reduction space and a basis defining this space.The functions β i , i=1, . . ., K, are called effective dimension reduction (edr)-directions and the space they generate is the edr-space.
FSIR is able to work with this model, indeed it yields d directions βi and corresponding (eigen-)values λi which allow one to rank the importance of βi , where βi and λi denote the estimated values.Obviously the β i (ν) are determined in the values ν j , j =1, . . ., d.The technique takes advantage of the method of inverse regression.Here the aim is not to estimate E[Y |R=r] but the reverse E[R|Y =y], a one-dimensional regression problem that avoids the curse of dimensionality.This curse means that high-dimensional spaces have too few data for local averaging, that is the risk or expected squared error of estimation of a nonparametric regression estimator increases rapidly with the dimension p and to maintain a given degree of accuracy of an estimator, the sample size must increase exponentially with the dimension p.This does not happen with parametric techniques, i.e. least squares regression, whose risk will decay to zero at a rate of n −1 , where n is the number of observations.With inverse regression instead of having one d-dimensional regression problem we have d one-dimensional regression problems which do not suffer from that curse.The key of FSIR is the connection between the edr-space and the inverse regression curve given by the covariance operator of the inverse regression curve e .Indeed it is possible to prove under some mild assumptions (Amato et al., 2006) that the eigenvalue-eigenvector decomposition of the operator −1 R e permits to identify a basis for the edr-space.Unfortunately the inverse of R is not bounded, therefore we consider ; in particular, we use the fact that has the same eigenvectors as R , where + e is the Moore-Penrose generalized inverse operator (Groetsch, 1974).
Let us now consider a dataset of N radiance spectra measured at wavenumbers ν j , j =1, . . ., d, and corresponding profiles Y N , that is (R n , Y n ), n=1, . . ., N , with R n being vectors of dimension d×1.Then after centering the data, R can be estimated by (2) For the estimate of e we may proceed in the following way.Let M Y (ν)=E(R|Y ) and MY (ν) be the wavelet smoothing of R(ν) with design points the Y n 's, n=1, . . ., N , obtained through the BINWAV estimator (Antoniadis and Pham, 1998).Then we consider the following estimate for e : with Convergence in probability of both ˆ R,N and ˆ e,N to R and e can be found in (Amato et al., 2006).To estimate them accurately, we improve the conditioning of ˆ e,N applying a projector method before performing the spectral decomposition: let π k N denote the orthogonal projector into the space spanned by the k N eigenvectors of ˆ e,N corresponding to the k N largest eigenvalues; we let ˆ k N e,N = πk N ˆ e,N π k N .Estimation of the EDR space is derived from the spectral decomposition of where is the pseudoinverse matrix defined through the Singular Value Decomposition (SVD) (Golub and Van Loan, 1983;Golub, 1970).Let (α i ) i=1,...,K be the smallest eigenvalues of (4) and η i the corresponding eigenfunctions, then Summarizing the procedure for computing an esti- and estimate the EDR directions by Eq. ( 5).
FSIR can be seen as a generalization of PCA.Indeed supposing that R is the identity operator, than FSIR aims at determining the directions along which to project the data by the eigenvalues-eigenvectors decomposition of e , which takes into account the information about the profiles by means of a regression on the spectra; on the contrary the PCA uses only spectral information.

Vertical resolution of the retrieval
If the vector, Ŷ =( Ŷ (1), Ŷ (2), . . ., Ŷ (M)) t , is made up with the layer-mean estimated values of a given parameter, in order to form the profile function of the parameter, then the retrieval covariance matrix can be obtained by where expectation value has to be taken with respect to training data set and Y true is the true value of the parameter.This matrix will be denoted by T , H 2 O and O 3 for temperature, water vapour and ozone profiles, respectively.The retrieval covariance matrix, for a given parameter profile, can be used to analyze the spatial vertical resolution of the parameter itself, according to a methodology which has been developed by Serio et al. (2008b) and which is summarized for the sake of clarity.
The vector (Y (1), Y (2), . . ., Y (M)) t represents the discretized version of a spatial function (i.e., temperature profile, water vapour mixing ratio profile, ozone mixing ratio U. Amato et al.: Dimension-reduction for IASI radiances profile).A strong correlation, that is a relatively high value of the covariance, between any two of the parameters means that the two have not been independently resolved by the data set, and that only some linear combination of the parameters is resolved.However, a direct examination and interpretation of the off-diagonal elements (covariances) of a covariance operator is not easy, which makes the definition of some suitable scalar index highly desirable.
To this end, it has to be considered that the retrieval correlation is determined by the non-null off-diagonal terms in the covariance operator.For a full independent retrieval (and therefore for a retrieval, which attains the maximum possible vertical spatial resolution), the covariance matrix has to be fully diagonal.
In what follows we will assume that the covariance operator has been normalized in order to obtain the correlation matrix, Ŷ may be additively decomposed in its diagonal and offdiagonal components: An index which assesses the dominance of the diagonal term over the off-diagonal one might be simply defined from the norms of the matrices in Eq. ( 8).However, the norm is not additive.In general, we have therefore, an index such as the ratio of the Ŷ , diag-norm to the Ŷ -norm is not well defined.It could be less or greater than one depending on the given matrix.
To quantify the relative contribution of the diagonal term (and off-diagonal term) to the norm of Ŷ , let us consider the SVD of Ŷ .With the usual notation, we have where U and V are unit, orthogonal matrices of size M by M; S is a diagonal matrix whose elements are the eigenvalues, positive definite and supposed decreasingly ordered, of the operator Ŷ .From the equation above, we have which with the position gives Finally, because of the definition of norm (e.g.Golub and Van Loan, 1983), we have where X(1, 1) is the element (1,1) of the given matrix X.The formula above is fully additive and allows us to decompose the norm of the retrieval covariance matrix in its diagonal and off-diagonal contribution.Thus, a proper index which quantifies the degree of diagonalization of Ŷ , that is how much the matrix is dominated by its diagonal terms, can be defined by For a full diagonal matrix, we have i D =1 and the retrieval is truly independent, whereas for a highly correlated matrix, we have i D =M −1 and we can retrieve only the columnar amount of a geophysical parameter.As far as IASI is concerned, this is, e.g., the case for trace gases such as carbon monoxide or methane.
The index (15) can be easily re-scaled to the range 1 to M by simply redefining it as Then, i D =1 simply means that for the retrieval at hand it is as if the full atmosphere had been divided just in one layer, that is only the columnar amount of the parameter has been resolved.On the opposite edge of the i D scale, we have i D =M, and the retrieval has been fully resolved on the grid mesh used to divide the atmosphere.Nearby layers can then, e.g., be used to form average quantities and, therefore, reduce the estimation error.

Application to IASI data
The retrieval accuracy of FSIR procedure has been assessed both on the Chevallier dataset and directly on IASI data.Analysis has focused on clear sky, sea-surface, tropical soundings.To limit the burden of the computational effort, only nadir view soundings have been considered.

Chevallier dataset
The training dataset consists of 377 (T , q, o) tropical profiles for clear-sky and sea surface extracted from the ECMWF compilation by Chevallier (Chevallier, 2001).
These profiles are the input to the σ -IASI code (Amato et al., 2002), which has been used to calculate the related synthetic IASI spectra at nadir view angle.
The spectral ranges used to train the regression algorithm include the intense CO 2 band 645 cm −1 -830 cm −1 , the ozone band 1010 cm −1 -1070 cm −1 , the window 1130 cm −1 -1180 cm −1 , a portion of the ν 2 water vapour band, 1400 cm −1 -1700 cm −1 , and the C 2 O/N 2 O shortwave absorption band, 2000 cm −1 -2230 cm −1 , for a total of 3305 channels.Thus, the data space is made up of radiance vectors of size d=3305.
The parameter space is made up of triplets (T j , q j , o j ) for each given atmospheric layer.The number of layers, M L , is the same as in σ -IASI, namely M L =60.For temperature we have one more parameter, since we consider the surface temperature, as well.
The FSIR is implemented according to a single parameter regression algorithm, in which each single element of the triplet (T j , q j , o j ) is regressed vs. the EDR components.
Similarly to the EOF methodology, the FSIR approach also requires that the user specifies the number of EDR components, p, (also called scores).To this purpose we define the root mean square error or estimation error, e(p), of the regression as e(p)=E[( Ŷ −Y Chev ) 2 ] 1/2 , where, as before, E[•] means expectation value, again Y denotes a generic parameter and Y Chev refers to values from the Chevallier dataset.Our algorithm is developed in the general context of a signal plus noise model (see, e.g., Eq. 1), therefore the curve e(p) reaches a sort of plateau or noise floor.A suitable choice for p is the closest value to the knee of the curve, e(p) (see Serio et al., 2009).Examples of these curves are shown in Fig. 1, which also provides a comparison with a PCA/EOF regression approach.Based on these curves the number of scores for which e(p) reaches a plateau ranges between p=5 and p=20.Also note that FSIR tends to be much more parsimonious in terms of the number of selected scores.
Figures 2-4 exemplify the expected root mean square error for temperature, water vapour and ozone, respectively.The figures also allow us to compare the forecast skill of FSIR with that of the EOF regression.For the temperature the two regression schemes are almost equivalent in terms of expected root mean square error.For water vapour FSIR is superior to EOF, the same as for ozone.

i D index
Using the methodology outlined in Sect.3, we can analyze the degree of interdependency of the retrieval for temperature, water vapour and ozone, indeed the index i D says how many independent linear combinations can be resolved from the data, and therefore how many independent pieces of information we have in the retrievals.From a statistical point of view it gives a measure of the degrees of freedom.
We give to the term degrees of freedom the usual meaning which is given to it in statistics.In its discretized form a given geophysical parameter is represented by a vector, v of size, say m.Therefore the vector v lives in a space of dimension m and we need, at least m independent equations, hence m pieces of information, to fully resolve for the geophysical parameter.However, because of data correlation and uncertainties, we have that the retrieved elements, vi can be correlated each to other, which means that only some linear combination of the elements of v has been resolved from the data set.The index i D says how many independent linear combina- tions can be resolved from the data, and therefore how many independent pieces of information we have in the retrieved, v.
For temperature the value of i D is 6.8 for FSIR and 5.6 for EOF.This means that FSIR gains more than one degree of freedom over the EOF regression.For water vapour and ozone we have 4.7 and 3.7, respectively, in the case of FSIR and 4.40 and 2.53 for EOF, respectively.The values for the index i D were evaluated on the whole training dataset consisting of 377 tropical profiles.
It is interesting to note that the FSIR retrieval is less correlated than that produced by the EOF scheme, which means that FSIR has a better capability to reveal features and structures along the vertical.However, in general the above i D values are quite small.For ozone, they indicate that only two or three pieces of information are available for the retrieval.For temperature and water vapour, i D says that only the very coarse features of the profile can be resolved.This is not a serious shortcoming for temperature, which is typically a smooth function of the altitude, but could become a serious limitations for water vapour, whose profile may be characterized by small-scale vertical structures.
In conclusion, FSIR can resolve important small structures in water vapour better than EOF.Of course, this does not mean that FSIR is, in absolute, the best retrieval method, we can only says that it is superior over EOF.

IASI retrieval and comparison with ECMWF analysis
The FSIR approach has been tested using IASI data obtained during the IASI commissioning phase.The cloud detection scheme described in Grieco et al. (2007) was applied to tropical spectra measured over the sea surface during the 22 July 2007.This yielded a total of 603 clear sky IASI spectra.To simplify the illustration of the results only nadir view soundings have been considered.
To develop a consistent set of truth data against which IASI retrieval could be compared, ECMWF atmospheric analysis fields for temperature, water vapour and ozone were considered.These fields where co-located in space and time to the 603 IASI soundings.We used atmospheric analysis fields at 00:00, 06:00, 12:00 and 18:00 UTC on 22 July 2007.At that time, the ECMWF model was characterized by a vertical discretization of the atmosphere into 60 pressure levels and a horizontal truncation of T511.This truncation corresponds to a grid spacing of about 40 km or, equivalently, to a horizontal grid box of 0.351 • ×0.351 • .The model has a hybrid vertical coordinate, with terrain-following coordinates in the lower troposphere and pressure coordinates in the stratosphere above about 70 hPa.Of the 60 levels in the vertical, 25 are above 100 hPa and the model top is at 0.1 hPa, corresponding to about 65 km.The vertical resolution of the analysis fields gradually decreases from 20 m at the surface to about 250 m at 1 km altitude, and about 1 km to 3 km in the stratosphere.The analysis fields were extracted from the ECMWF archive at the full T511 resolution, interpolated to a grid of points with a separation of 0.3 • ×0.3 • and then colocated to the IASI soundings.The statistics of the difference between global radiosonde observations and ECMWF analysis in the troposphere show values of the standard deviation typically between 0.5 and 1 K for temperature and between 0.5 and 1.5 g/kg for water vapour.In addition to fields of temperature, water vapour and ozone, ECMWF fields of seasurface temperature (SST) were also used in the study.It should be noted that these fields are based on analyses received daily from the National Centers for Environmental Prediction (NCEP), Washington DC, on a 0.5 • ×0.5 • grid.These analyses are based on ship, buoy and satellite observations.In shallow waters, where rapid changes due to the upwelling radiation can occur close to land, the observed SST can sometimes differ as much as 5 K from the NCEP analysis.
One important aspect which has been possible to analyze with the help of the ECMWF analysis is the sensitivity of the retrieval accuracy to the choice of the number, p, of FSIR and/or PCA scores.Because of many factors, which include radiative transfer accuracy, noise specifications, cloud contamination, in practice it may happen that the optimal choice, p opt , performed based on the training data set may result to be sub-optimal.Thus, when this value is applied to real data and real conditions, a misfit may occur, which is much larger than that expected on the basis of the training data set.
To address this point we have redefined the root mean square difference, e(p), as e(p)=E[( Ŷ −Y ecmwf ) 2 ] 1/2 .For temperature, e(p) is shown in Fig. 5 for three values of p, namely p opt , p opt +5 and p opt +10.The calculations have  been performed for FSIR and PCA.It can be seen that, in comparison to EOF, FSIR is more robust in terms of accuracy to variations in the number of scores, p.This is much more evident for the case of water vapour (Fig. 6) and ozone (Fig. 7).Apart from some isolated error spikes, evident in the case of ozone, FSIR also provides a more accurate retrieval, once compared to EOF.This can be seen in the three next Figs.8-10 which compare FSIR performance to EOF for the case of temperature, water vapour and ozone.
Figures 8-10 also suggest that FSIR accuracy for temperature is within 1-2 K, which is a bit larger than the 1 K expected for IASI.However, it is important to consider that the ECMWF analysis also has an uncertainty, which, as discussed above, is of the order of 0.5 to 1 K. Thus in e(p) we also include the uncertainty of the ECMWF analysis.
The same as above can be said for water vapour: including the ECMWF uncertainty, the performance of FSIR reaches the expected accuracy of 10% only in the very deepest part of the atmosphere.
For ozone we get a very smooth retrieval, as it was expected based on the i D index for this parameter.However, the result compares fairly well with ECMWF mean profile.

Conclusions
A new statistical strategy based on dimension reduction for the retrieval of atmospheric parameters from IASI radiances has been presented and discussed.Applications to IASI data have been considered for the case of tropical soundings.The new strategy, FSIR, has been also compared to a usual EOF regression scheme.The comparison shows that, mostly for gases (we have analyzed water vapour and ozone), FSIR gets a higher performance with respect to EOF.In general, FSIR seems to provide a retrieval which is better resolved along the vertical, which is particularly interesting for water vapour.
For temperature the FSIR scheme provides results quite close to the expected performance of 1 K in the lower part of the atmosphere.For water vapour the goal of 10% accuracy is reached only for the very lower part of the atmosphere.
For ozone it seems that FSIR is capable to provide at least 3 pieces of information.
We think that FSIR can provide a valid initialization scheme for physical-based inversion strategy, which is now under investigation.

Fig. 1 .
Fig. 1.Root mean square error of the retrievals as a function of the number of scores, for various atmospheric layers and parameters.

Fig. 2 .
Fig. 2. Temperature root mean square error of the retrievals on the Chevallier dataset as obtained for FSIR and EOF.

Fig. 3 .Fig. 4 .
Fig. 3. Water vapour percentage root mean square error of the retrievals on the Chevallier dataset as obtained for FSIR and EOF.

Fig. 7 .
Fig. 7. Ozone root mean square difference (IASI retrieval -ECMWF) for 3 choices of the number of EDR (PC) scores for FSIR (left plot) and EOF (right plot).
Fig. 9. (a) -Mean retrieved water vapour profile obtained by averaging over the 603 IASI soundings and comparison with the ECMWF corresponding mean profile.(b) Percentage root mean square difference (IASI retrieval -ECMWF) as obtained for EOF and FSIR methodologies.
Fig. 10.(a) -Mean retrieved ozone profile obtained by averaging over the 603 IASI soundings and comparison with the ECMWF corresponding mean profile.(b) Percentage root mean square difference (IASI retrieval -ECMWF) as obtained for EOF and FSIR methodologies.