The SpatialExtremes Package: An R package to model spatial extremes

The SpatialExtremes package aims to provide tools for modeling spatial extremes. Many environmetal processes are spatial by nature and the knowledge of the univariate/multivariate distribution of extremes is not enough. This package aims to fill this gap by using several approaches. Currently, the modeling is performed through the max-stable framework. Max-stable processes are the extension of the extreme value theory to random fields. Later, different approaches will be implemented such as the latent variable and the copula based methodologies.
The development of the package has been financially supported by the Competence Center Environment and Sustainability (CCES) and more precisely within the EXTREMES project.

Features

Currently, the SpatialExtremes package models spatial extremes by fitting max-stable processes to data using composite likelihood maximisation. Many tools to characterize the spatial behaviour are also implemented.
In a nutshell, the SpatialExtremes packages can:

Fit max-stable processes to data
Estimate the extremal coefficient
Allow for prediction
Fits penalized splines with radial basis functions
Plot the covariance and the madogram to better understand how evolve dependence in space

The SpatialExtremes Package in a Few Lines

In this section, we explicit some of the most useful functions of the package. However, for a full description, users may want to have a look to the package vignette and the html help of the package.

Fit a max-stable process - assuming unit Frechet margins:

##Fit the the Smith's model to data:
fitmaxstab(data, coord, "gauss", fit.marge = FALSE) Estimator: MPLE
Model: Smith
Pair. Deviance: 1049771
TIC: 1049967
Covariance Family: Gaussian

Estimates
Marginal Parameters:
Not estimated.
Dependence Parameters:
cov11 cov12 cov22
10.175 5.305 4.236

Standard Error Type: score

Standard Errors
cov11 cov12 cov22
1.0323 0.5549 0.3170
...
##Fit the Schlather's model with the powered exponential covariance function:
fitmaxstab(data, corrd, "powexp", fit.marge = FALSE) Estimator: MPLE
Model: Schlather
Pair. Deviance: 154654.7
TIC: 154744
Covariance Family: Powered Exponential

Estimates
Marginal Parameters:
Assuming unit Frechet.
Dependence Parameters:
sill range smooth
0.790 17.519 1.664

Standard Error Type: score

Standard Errors
sill range smooth
0.02183 3.05162 0.24029
...

Fit a max-stable process - allowing for common GEV margins

##With linear models for the GEV parameters:


      fitmaxstab(data, coord, "gauss", loc.form = y ~ lon, scale.form
      = y ~ lat, shape.form = y ~ 1)

Estimator: MPLE
Model: Smith
Pair. Deviance: 325981.1
TIC: 16548523
Covariance Family: Gaussian

Estimates
Marginal Parameters:
Location Parameters:
locCoeff1 locCoeff2
5.560 1.962
Scale Parameters:
scaleCoeff1 scaleCoeff2
3.9727 0.9228
Shape Parameters:
shapeCoeff
0.2017
Dependence Parameters:
cov11 cov12 cov22
8.862 4.187 3.470

Standard Error Type: score
...
##With a p-spline on the location parameter:


      loc.form <- y ~ rb(lon, knots, degree, penalty)

      fitmaxstab(data, coord, "gauss", loc.form = loc.form, scale.form =
      y ~ lat, shape.form = y ~ 1)

Fit a p-spline

##Setting up the smoothing parameter using CV/GCV:


      cv(y, x, knots, degree); gcv(y, x, knots, degree);

$cv
[1] 91.0042

$penalty
[1] 0.1914471

$nlm.code
[1] 1
##Fit a p-spline using different smoothing parameters:


      rbpspline(y, x, knots, degree, penalty = "gcv");

      rbpspline(y, x, knots, degree, penalty = "cv");

      rbpspline(y, x, knots, degree, penalty = .2);

Call:
rbpspline(y = y, x = x, knots = knots, degree = 3)

Rank: 200 (G)CV Score: 0.43
Degree: 3 Penalty: 0.221

Degree of freedom: 9.856
Res. Degree of freedom: 187.49

Model selection

##Using the TIC as an analogue of AIC under misspecification:


      TIC(fitted1, fitted2)

fitted1 fitted2
448784.9 448839.5
##Anova (nested models) under misspecifications:


      anova(fitted1, fitted2); anova(fitted1, fitted2, method = "CB")

Eigenvalue(s):
183.4

Analysis of Variance Table
MDf Deviance Df Chisq Pr(> sum lambda Chisq)
M2 1 66982
M1 2 64155 1 2826.6 8.642e-05 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Plots

##Covariance functions:
covariance(sill = 1, range = 1, smooth = 0.5, cov.mod = "whitmat");

##Madogram:
madogram(data, coord, which = "mado")

##Extremal coefficient:
extcoeff(fitted.object);

##Prediction of return levels:

map(fitted.object, "quant", ret.per =
    50)

condmap(schlather, c(1, 1), seq(0, 10, length
    = 20), seq(0,10, length = 20))

##Semi-parametric estimates of the extremal coefficient:

fitextcoeff(data, coord); fitextcoeff(data, coord, angles =
    seq(-pi, pi, length = 4))

##Profile composite likelihood:
profile(fitted.object, "cov11", range = c(80, 200))

##Plot a p-spline:


      fitted <- rbpspline(y, x, knots, degree); lines(fitted, col =
      2)

Manuals

We have written a package vignette to help new users. This user's guide is a part of the package - just run vignette("SpatialExtremesGuide") in an R console once the package is loaded.

Contribute to the Project

If you are interested in joining the project, you must first create an R-forge account. Then, just go here.

Contact

Any suggestions, feature requests, bugs: select the appropriate tracker
Author: Mathieu Ribatet (homepage)