Franco Villoria, Maria

Ignaccolo, Rosaria

English

FRANCO VILLORIA, Maria

Institutional Research Information System University of Turin

This is an author version of the contribution published on:Rosaria Ignaccolo, Maria Franco VilloriaKriging for functional data: uncertainty assessmentEditor: CUEC Cooperativa Universitaria Editrice Cagliaritana2014ISBN: 9788884678744in47th Scientific Meeting of the Italian Statistical Society P R O C E E D I N GS1 - 647th SIS Scientific Meeting of the Italian Statistica SocietyCagliariJune 2014Kriging for functional data: uncertaintyassessmentKriging per dati funzionali: valutazione dell’incertezzaMaria Franco-Villoria and Rosaria IgnaccoloAbstract We predict a curve at an unmonitored site taking into account exogenousvariables using a functional kriging model with external drift and, alternatively, anadditive model with a spatio-temporal smooth term. To evaluate uncertainty of thepredicted curves, a semi-parametric bootstrap approach is used for the first, whilestandard inference is used for the second. The performance of both approaches isillustrated on pollutant functional data.Abstract Allo scopo di predire una curva in un sito non monitorato, tenendo contodi variabili esogene, viene usato un modello di kriging funzionale con deriva esternaed in alternativa un modello additivo con una componente spazio-temporale.L’incertezza viene valutata usando un approccio bootstrap semi-parametrico nelprimo caso e inferenza classica nel secondo. La performance dei due approcci e´illustrata su dati funzionali di inquinamento.Key words: kriging, P-spline, confidence bands, performance index1 IntroductionSpatial prediction of particulate matter (PM) concentration is useful in order to as-sess air quality and health risk where no monitoring stations are available. The ob-served air pollutants and meteorological variables time series, collected at variouslocations of a monitoring network, can be treated as spatially dependent functionaldata. Several papers consider ordinary kriging models for functional data under theassumption of a constant mean (e.g. [4, 5, 9]). The more recent ones consider theMaria Franco-Villoria, Rosaria IgnaccoloDipartimento di Economia e Statistica, Universita` degli Studi di Torino, Italy, e-mail:maria.francovilloria@unito.it, rosaria.ignaccolo@unito.itWork partially supported by FIRB 2012 grant (project no.RBFR12URQJ) provided by the ItalianMinistry of Education, Universities and Research.12 Maria Franco-Villoria and Rosaria Ignaccolomean as a function of longitude and latitude [1, 8] and (both scalar and functional)exogenous variables [6]. However, uncertainty evaluation about prediction remainsan open issue, as kriging variance is constant in time. Given the difficulty to de-rive sampling distributions for functional data, confidence band calculation can beapproached using resampling methods. We adapt the semi-parametric bootstrap ap-proach for spatially correlated data proposed in [7] to the case of functional data.Confidence bands are obtained by ordering the bootstrapped predicted curves in twodifferent ways, based on functional depth and on distance between curves. Alterna-tively, spatially dependent functional data can be modelled in a longitudinal per-spective by means of an additive model that includes a smooth function of longitudeand latitude, exploiting the close connection between penalized spline smoothingand kriging [12]. Uncertainty for a predicted curve can be assessed using classicalinference. The two approaches are illustrated on pollutant functional data.2 Case study: PM10 concentration in PiemonteOur case study consists of daily PM10 concentrations (in µg/m3) measured fromOctober 2005 to March 2006 by the monitoring network of Piemonte region (Italy)in 24 sites (red triangles in Fig. 1(a)), and in 10 more sites used as validation stations(blue dots in Fig. 1(a)). Data were log transformed to stabilize the within-stationvariances and normalize the marginal distribution of PM10 data. The covariates con-sidered are: coordinates and altitude (scalar); daily maximum mixing height, dailytotal precipitation, daily mean wind speed, daily mean temperature and daily emis-sion rates of primary aerosols (functional). Since the ranges of the covariates arequite different, a standardization procedure is applied (for further details see [2]).3 Functional Kriging with External Drift (FKED)Letϒs = {Ys(t);t ∈ T} be a functional random variable observed at location s∈D⊆Rd , whose realization is a function of t ∈ T , where T is a compact subset of R.Assume that we observe a sample of curvesϒsi , for si ∈D, i = 1, . . . ,n, that take val-ues in a separable Hilbert space of square integrable functions. The set {ϒs,s ∈ D}constitutes a functional random field or a spatial functional process [4], that canbe non-stationary and that is supposed to be decomposed as ϒs = µs + εs. The termµs is interpreted as a drift describing a spatial trend while εs represents a residualrandom field that is zero-mean, second-order stationary and isotropic. At site si,i = 1, . . . ,n, and point t, the model can be rewritten as a functional concurrent linearmodel Ysi(t) = µsi(t)+ εsi(t) with the driftµsi(t) = α(t)+∑pγp(t)Cp,i +∑qβq(t)Xq,i(t) (1)Kriging for functional data: uncertainty assessment 3where α(t) is a functional intercept, Cp,i is the p-th scalar covariate at site si, Xq,i isthe q-th functional covariate at site si, γp(t) and βq(t) are the covariate coefficientsand εsi(t) represents the residual spatial functional process {εs(t), t ∈ T,s ∈ D} atsite si. Model (1) is fitted by means of a GAM representation (see [6]), and thefunctional residuals esi(t) = Ysi(t)− µˆsi(t) can be used to predict the residual curveat a new site s0 via ordinary kriging for functional data [5], according to whicheˆs0(t) = ∑ni=1 ξiesi(t), with kriging coefficients ξi ∈ R. More complex alternativesconsidering non-constant kriging coefficients are available [6]. Prediction at the un-monitored site s0 is obtained by adding up the two terms, i.e. Yˆs0(t) = µˆs0(t)+ eˆs0(t),where µˆs0(t) = αˆ(t)+∑p γˆp(t)Cp,0+∑q βˆq(t)Xq,0(t).3.1 Uncertainty Evaluation: a bootstrap approachTo evaluate the uncertainty of a predicted curve Yˆs0(t), t = 1, . . . ,M, we extend thesemi-parametric bootstrap approach for spatially correlated data proposed by [7] tothe functional context. The bootstrapping algorithm can be summarized as follows:1. Estimate and remove the drift µs following Model (1) to obtain the functionalresiduals esi(t) = Ysi(t)− µˆsi(t).2. Estimate the residuals covariance matrix Σ through the trace-semivariogram:υˆ(h) =12 |N(h)| ∑i, j∈N(h)∫T(esi(t)− es j(t))2dtwhere N(h) = {(si,s j) :∥∥si− s j∥∥ = h}. A parametric model (exponential orMate´rn for example) can be fitted to the points (hg, υˆ(hg)),g = 1, . . . ,G, as itis done in classical geostatistics. Using the Cholesky decomposition Σˆ = UˆUˆT ,the functional residuals can be transformed to become spatially uncorrelated:ζn×M =(ζ (s1)′, . . . ,ζ (sn)′)′= Uˆ−1n×n (Yn×M − µˆn×M) .3. Generate B bootstrap samples ζ ∗(s1), . . . ,ζ∗(sn) from ζ (s1), . . . ,ζ (sn) using thesmoothed bootstrap as suggested in [3], replacing the empirical distribution func-tion of {ζ (s1), . . . ,ζ (sn)}, denoted as Fn, by a smooth version Fˆn .4. The final bootstrap sample is determined using an inverse transform:Y ∗si (t) = µˆsi(t)+Uˆζ∗si(t).The bootstrap samples are then fed into the FKED method to obtain B predictioncurves at the 10 validation sites. These B prediction curves need to be ordered todetermine the upper and lower limits of the confidence band. Here we consider twoordering techniques based on band depth and L2 distance between curves.Band depth [11] can be defined for any set of k curves (here k = 2). The sampleband depth (BD) of y(t) can be calculated as the proportion of bands delimited by 24 Maria Franco-Villoria and Rosaria Ignaccolocurves containing the whole curve y(t) [11]. The modified band depth (MBD), thattakes into account whether a portion of the curve is in the band, is defined asMBDn,2(y)=(n2)−1∑1≤i1<i2≤nλ ({t ∈ T : minr=i1,i2yr(t)≤ y(t)≤ maxr=i1,i2yr(t)})λ (T )where λ is the Lebesgue measure on T (for details see [11]). The lower/upper limitsof a 95% confidence band (based on band depth) are obtained by taking the point-wise (w.r.t. t) minimum/maximum of the 95% deepest curves (i.e. closest to thecenter of the distribution). On the other hand, the 95% confidence ball (based on L2distance) is made of the 95% curves closest to the FKED predicted curve Yˆs0(t) [3].3.2 Uncertainty Evaluation: a GAM approachWe incorporate a smooth function of longitude, latitude and time, namely f (lon, lat, t),in Model (1). By setting a penalized bivariate spline basis for longitude and latitudea spatial covariance structure is implicit in the model [10], allowing for spatial pre-diction. The model is linear in the coefficients and can be written in matrix formas Yˆ = SY where S = X(X ′X +ηP)−1X ′ is the smoothing matrix, X is the designmatrix, P is the penalty matrix and η is the smoothing parameter. At a new locations0, the predicted value is given by Yˆs0(t) = Ss0y, where Ss0 = Xs0(X′X +ηP)−1X ′.Approximate 95% confidence bands can be calculated as [12]:Yˆs0(t)± 1.96σˆε√1+ ||Ss0||2. (2)Alternatively, uncertainty can be evaluated using the reformulation of the model asa mixed model so that bias is taken into account (see [12], p.138).4 ResultsPredictions at the 10 validation sites were obtained by means of FKED as summa-rized in Section 3 and a fully additive model as described in Section 3.2. Four perfor-mance indices were calculated to compare the spatial prediction of the two modelsconsidered and are summarized in Figure 1(b). The indices are: RMSE, Pearsoncorrelation, Normalized Mean Bias Factor and weighted normalized MSE of thenormalized ratios (WNNR), that takes into account peaks in observed data, wherethe normalized ratio is defined as exp(−| ln( f itted/observed)|) (see [6] for details).All four indices suggest that FKED provides better predictions for the 10 validationsites. For each of these sites, a bootstrap sample of predicted curves of size 1000was obtained following the algorithm illustrated in Section 3.1. Band depth wascalculated using the modified version MBD. The resulting confidence balls/bands,Kriging for functional data: uncertainty assessment 5as well as 95% confidence bands according to Eq. (2), are shown in Figure 2 fortwo of the 10 sites, namely 25-Biella-Largo Lamarmora and 30-Saliceto. Overall,the two (band and ball) bootstrap based uncertainty measures seem to agree well.For Biella-Largo Lamarmora, nearly all (92.3%) the observed data lie within thebootstrap based confidence ball/band, while for Saliceto, there is a high percentage(32.4%) of observed values that lie outside the confidence ball/band. This reflectsthe performance of the corresponding kriging predictions, good in the first valida-tion site (the predicted curve agrees well with the observed values), but not in thesecond. On the other hand, all data points lie outside the inference confidence bands,which are much narrower than the bootstrap ones. Part of this difference can be at-tributed to the fact that bootstrap based confidence bands are simultaneous, whileinference based ones are pointwise.5 Discussion and Future WorkSpatial functional analysis provides an alternative to spatio-temporal modelling ca-pable of predicting a whole curve taking into account exogenous covariates. How-ever, uncertainty evaluation remains an open issue. Here, two different approachesare illustrated. The resulting confidence bands differ considerably, with the boot-strap based ones being too wide while inference based ones are too narrow. Despitebeing easy to implement, the bootstrap approach seems to be non ideal; it can becomputationally expensive and in our case study there are time periods at whichthe confidence region becomes unjustifiably wider. In the case of the fully additivemodel, the width of the (inference based) interval remains fairly constant, beingslightly wider at the beginning and end of the predicted curve. Improving uncer-tainty estimation through a mixed model approach is part of our ongoing research.References1. Caballero, W., Giraldo, R., Mateu, J.: A universal kriging approach for spatial functional data.Stochastic Environmental Research and Risk Assessment 27(7), 1553–1563 (2013)2. Cameletti, M., Ignaccolo, R., Bande, S.: Comparing spatiotemporal models for particulatematter in Piemonte. Environmetrics 22, 985-996 (2011)3. Cuevas, A., Febrero, M., Fraiman, R.: On the use of the bootstrap for estimating functionswith functional data. Comput. Statist. Data Anal. 51, 1063–1074 (2006)4. Delicado, P., Giraldo, R., Comas, C., Mateu, J.: Statistics for spatial functional data: somerecent contributions. Environmentrics 21, 224–239 (2010)5. Giraldo, R., Delicado, P., Mateu, J.: Ordinary kriging for function-valued spatial data. Envi-ronmental and Ecological Statistics 18(3), 411–426 (2011)6. Ignaccolo, R., Mateu, J., Giraldo, R.: Kriging with external drift for functional data for airquality monitoring. Stochastic Environmental Research and Risk Assessment (2013) doi:10.1007/s00477-013-0806-y6 Maria Franco-Villoria and Rosaria Ignaccolo(a)FKED GAM0.20.30.40.50.60.7RMSEFKED GAM0.600.650.700.750.800.850.900.95PearsonFKED GAM−0.10−0.050.000.05NMBFFKED GAM0.010.020.030.04WNNR(b)Fig. 1 (a) Locations of 24 PM10 monitoring sites (red triangles) and 10 validation stations (bluedots). (b) Boxplot of performance indices distribution over the validation stations.Fig. 2 Original data (blackdots), GAM predicted curve(green) and corresponding95% confidence band (dashedgreen line), FKED predictedcurve (red), 95% confidenceball (grey) based on L2 dis-tance and 95% confidenceband (blue) based on MBDfor locations 25 and 30 inFigure 1(a).0 50 100 1501234525 − Biella Largo LamarmoraDaylog(PM10)0 50 100 1501234530 − SalicetoDaylog(PM10)7. Iranpanah, N., Mohammadzadeh, M., Taylor, C.C.: A comparison of block and semi-parametric bootstrap methods for variance estimation in spatial statistics. Comput. Statist.Data Anal. 55, 578–587 (2011)8. Menafoglio, A., Secchi, P., Dalla Rosa, M.: A Universal Kriging predictor for spatially depen-dent functional data of a Hilbert Space. Electronic Journal of Statistics 7, 2209–2240 (2013)9. Nerini, D., Monestiez, P., Mant, C.: Cokriging for spatial functional data. J. Multivariate Anal.101, 409-418 (2010)10. Nychka, D.: Spatial process estimates as smoothers. In: Schimek, M.G. (ed.), Smoothingand Regression. Approaches, Computation and Application, pp. 393-424. Wiley, New York(2000)11. Lopez-Pintado, S., Romo, J.: On the concept of depth for Functional Data. J. Amer. Statist.Assoc. 104(486), 718–734 (2009)12. Ruppert, D., Wand, M.P., Carroll, R.J.: Semiparametric Regression. Cambridge UniversityPress, New York (2003)

Kriging for functional data: uncertainty assessment

https://iris.unito.it/bitstream/2318/145817/1/2014_SIS_fked_2953_4aperto_1343014.pdf

Kriging for functional data: uncertainty assessment

Abstract

Similar works

Full text

Available Versions

Institutional Research Information System University of Turin