Title: | Blind Source Separation for Multivariate Spatio-Temporal Data |
---|---|
Description: | Simultaneous/joint diagonalization of local autocovariance matrices to estimate spatio-temporally uncorrelated random fields. |
Authors: | Christoph Muehlmann [aut] , Nikolaus Piccolotto [aut] , Claudia Cappello [aut] , Sandra De Iaco [aut] , Klaus Nordhausen [aut, cre] |
Maintainer: | Klaus Nordhausen <[email protected]> |
License: | GPL (>= 2) |
Version: | 0.4-0 |
Built: | 2024-11-27 06:29:26 UTC |
Source: | CRAN |
Simultaneous/joint diagonalization of local autocovariance matrices to estimate spatio-temporally uncorrelated random fields.
Package: | SpaceTimeBSS |
Type: | Package |
Version: | 0.4-0 |
Date: | 2024-05-29 |
License: | GPL (>= 2) |
Solving the second order blind source separation problem for multivariate space-time random fields. The random fields can be irregular in space but must be regular in time.
The main function of this package is:
stbss
This function computes local autocovariance matrices. The considered temporal lags are integer numbers and the spatial lags are defined by spatial kernel functions. Then, these local autocovariance matrices and the sample covariance are simultaneously/jointly diagonalized.
The package also contains a 9-variate dataset of deseasonalized weekly climate and meteorological measurements from the Italian Veneto region between 2000 and 2022 meteo_veneto
.
Joint diagonalization is computed with the frjd
(fast real joint diagonalization) algorithm from the package JADE
.
The available finite realizations of the space time random fields can be defined by matrices or an object of classes STFDF
, STSDF
or st_sftime
.
Christoph Muehlmann, Nikolaus Piccolotto, Claudia Cappello, Sandra De Iaco, Klaus Nordhausen
Maintainer: Klaus Nordhausen [email protected]
Muehlmann, C., De Iaco, S. and Nordhausen, K. (2023), Blind Recovery of Sources for Multivariate Space-Time Environmental Data. Stochastic and Environmental Research and Risk Assessment, 37, 1593–1613, <doi:10.1007/s00477-022-02348-2>.
Extracts the estimated unmixing matrix of an object of class 'stbss'
.
## S3 method for class 'stbss' coef(object, ...)
## S3 method for class 'stbss' coef(object, ...)
object |
object of class |
... |
further arguments to be passed to or from methods. |
Returns the estimated unmixing matrix of an object of class 'stbss'
as a numeric matrix. For a description of the class 'stbss'
see stbss
.
Computation of local autocovariance matrices for a multivariate space-time dataset based on a given set of spatio-temporal kernel functions.
lacov(x, coords, time, kernel_type, kernel_parameters, lags, kernel_list = NULL, center = TRUE)
lacov(x, coords, time, kernel_type, kernel_parameters, lags, kernel_list = NULL, center = TRUE)
x |
either a numeric matrix of dimension |
coords |
a numeric matrix of dimension |
time |
a numeric vector of length |
kernel_type |
either a string or a string vector of length |
kernel_parameters |
a numeric vector of length |
lags |
an integer vector of length |
kernel_list |
a list of spatio-temporal kernel matrices with dimension |
center |
logical. If |
Local autocovariance matrices are defined by
with
Here, are the
p
random field values at location ,
is the sample mean vector, and the space-time kernel function
determines the locality. The following kernel functions are implemented and chosen with the argument
kernel_type
:
'ring'
: the spatial parameters are inner radius and outer radius
, with
, and
, the temporal parameter is the temporal lag
:
'ball'
: the spatial parameter is the radius , with
, the temporal parameter is the temporal lag
:
'gauss'
: Gaussian function where 95% of the mass is inside the spatial parameter , with
, the temporal parameter is the temporal lag
:
Above, represents the indicator function. The argument
kernel_type
determines the used kernel function as presented above, the argument lags
provides the used temporal lags for the kernel functions ( in the above formulas) and the argument
kernel_parameters
gives the spatial parameters for the kernel function. Each of the arguments kernel_type
, lags
and kernel_parameters
can be of length K
or 1
. Specifically, kernel_type
can be either one kernel, then each local autocovariance matrix use the same kernel type, or of length K
which leads to different kernel functions for the provided kernel parameters. lags
can be either one integer, then for each kernel the same temporal lag is used, or an integer vector of length K
which leads to different temporal lags. In the same fashion kernel_parameters
is a vector of length K
or 1
. If kernel_type
equals 'ball'
or 'gauss'
then the corresponding entry of kernel_parameters
gives the single spatial radius parameter. In contrast, if (at least one entry of) kernel_type
equals 'ring'
, then kernel_parameters
must be a list of length K
(or 1
) where each entry is a numeric vector of length 2
defining the inner and outer spatial radius. See examples below.
Alternatively, a list of kernel matrices can be given directly to the function lacov
through the kernel_list
argument. A list with kernel matrices can be computed with the function stkmat
.
lacov
returns a list of length K
where each entry is a numeric matrix of dimension c(p, p)
corresponding to a local autocovariance matrix.
Muehlmann, C., De Iaco, S. and Nordhausen, K. (2023), Blind Recovery of Sources for Multivariate Space-Time Environmental Data. Stochastic and Environmental Research and Risk Assessment, 37, 1593–1613, <doi:10.1007/s00477-022-02348-2>.
# space and time coordinates n_t <- 50 n_sp <- 10 st_coords <- as.matrix(expand.grid(1:n_sp, 1:n_sp, 1:n_t)) # simulate three latent white noise fields field_1 <- rnorm(nrow(st_coords)) field_2 <- rnorm(nrow(st_coords)) field_3 <- rnorm(nrow(st_coords)) # compute the observed field latent_field <- cbind(field_1, field_2, field_3) mixing_matrix <- matrix(rnorm(9), 3, 3) observed_field <- latent_field # lacov with different ring kernels and same lags lacov_r <- lacov(observed_field, coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ring', kernel_parameters = list(c(0, 1), c(1, 2)), lags = 1) # lacov with same ball kernels and different lags lacov_b <- lacov(observed_field, coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ball', kernel_parameters = 1, lags = c(1, 2, 3)) # lacov with different gauss kernels and different lags lacov_g <- lacov(observed_field, coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'gauss', kernel_parameters = 1, lags = 1:3) # lacov mixed kernels lacov_m <- lacov(observed_field, coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = c('ball', 'ring', 'gauss'), kernel_parameters = list(1, c(1:2), 3), lags = 1:3) # lacov with a kernel list kernel_list <- stkmat(coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ring', kernel_parameters = list(c(0, 1)), lags = 1) lacov_k <- lacov(observed_field, kernel_list = kernel_list)
# space and time coordinates n_t <- 50 n_sp <- 10 st_coords <- as.matrix(expand.grid(1:n_sp, 1:n_sp, 1:n_t)) # simulate three latent white noise fields field_1 <- rnorm(nrow(st_coords)) field_2 <- rnorm(nrow(st_coords)) field_3 <- rnorm(nrow(st_coords)) # compute the observed field latent_field <- cbind(field_1, field_2, field_3) mixing_matrix <- matrix(rnorm(9), 3, 3) observed_field <- latent_field # lacov with different ring kernels and same lags lacov_r <- lacov(observed_field, coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ring', kernel_parameters = list(c(0, 1), c(1, 2)), lags = 1) # lacov with same ball kernels and different lags lacov_b <- lacov(observed_field, coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ball', kernel_parameters = 1, lags = c(1, 2, 3)) # lacov with different gauss kernels and different lags lacov_g <- lacov(observed_field, coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'gauss', kernel_parameters = 1, lags = 1:3) # lacov mixed kernels lacov_m <- lacov(observed_field, coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = c('ball', 'ring', 'gauss'), kernel_parameters = list(1, c(1:2), 3), lags = 1:3) # lacov with a kernel list kernel_list <- stkmat(coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ring', kernel_parameters = list(c(0, 1)), lags = 1) lacov_k <- lacov(observed_field, kernel_list = kernel_list)
Weekly aggregated climate and meteorological deseasonalized data in Veneto region (Italy) for a 23-year span (2000-2022).
data("meteo_veneto")
data("meteo_veneto")
Object of class data.frame with 85248 rows and 13 variables, where 85248 consists of 1184 weekly observations times 72 spatial locations. The variables are as follows:
x
x coordinates in meters (Gauss Boaga - EPSG:3003)
y
y coordinates in meters (Gauss Boaga - EPSG:3003)
sp.ID
code for the spatial location
timeIndex
code for the temporal observation, from 1 to 1184
deseas_ET0
evapotranspiration levels (mm)
deseas_rad
solar radiation (MJ/m^2)
deseas_tmax
maximum temperature (degrees C)
deseas_taver
average temperature (degrees C)
deseas_tmin
minimum temperature (degrees C)
deseas_hmax
maximum humidity (%)
deseas_hmin
minimum humidity (%)
deseas_wind
wind velocity (m/s)
deseas_log_prec
log of precipitation (mm)
The evapotranspiration levels were estimated by ARPA Veneto according to the Hargreaves equation. The data have been obtained by removing the annual periodicity from the raw data and then computing weekly averages.
The raw data can be downloaded from the Environmental Protection Agency of Veneto Region (ARPA Veneto) website.
Prints the estimated unmixing matrix, the (pseudo-)eigenvalues and the diagonalized local autocovariance matrices for an object of class 'stbss'
.
## S3 method for class 'stbss' print(x, ...)
## S3 method for class 'stbss' print(x, ...)
x |
object of class |
... |
additional arguments for the method |
No return value.
For a given multivariate space-time dataset, stbss
estimates the realization of spatio-temporally uncorrelated random fields through a linear transformation which is defined by a so-called mixing matrix and a location vector. This is done assuming a spatio-temporal blind source separation model and simultaneously/jointly diagonalizing the sample covariance matrix and one/many local autocovariance matrices.
stbss(x, ...) ## Default S3 method: stbss(x, coords, time, kernel_type, kernel_parameters, lags, ordered = TRUE, kernel_list = NULL, ...) ## S3 method for class 'STFDF' stbss(x, ...) ## S3 method for class 'STSDF' stbss(x, ...) ## S3 method for class 'sftime' stbss(x, ...)
stbss(x, ...) ## Default S3 method: stbss(x, coords, time, kernel_type, kernel_parameters, lags, ordered = TRUE, kernel_list = NULL, ...) ## S3 method for class 'STFDF' stbss(x, ...) ## S3 method for class 'STSDF' stbss(x, ...) ## S3 method for class 'sftime' stbss(x, ...)
x |
either a numeric matrix of dimension |
coords |
a numeric matrix of dimension |
time |
a numeric vector of length |
kernel_type |
either a string or a string vector of length |
kernel_parameters |
a numeric vector of length |
lags |
an integer vector of length |
ordered |
logical. If |
kernel_list |
a list of spatio-temporal kernel matrices with dimension |
... |
further arguments for the fast real joint diagonalization algorithm that jointly diagonalizes the local covariance matrices. See details and |
It is assumed that the p-variate space-time random field is formed by
where is the latent p-variate space-time random field,
and
are the mixing matrix and a location vector and
and
are the space and time coordinates. Furthermore, it is assumed that
is white and consists of space-time uncorrelated components. The goal is to reverse the linear form by estimating an unmixing matrix and the location vector. This is done by simultaneously/jointly diagonalizing local autocovariance matrices which are defined by
with
Here, are the
p
random field values at location ,
is the sample mean vector, and the space-time kernel function
determines the locality. The following kernel functions are implemented and chosen with the argument
kernel_type
:
'ring'
: the spatial parameters are inner radius and outer radius
, with
, and
, the temporal parameter is the temporal lag
:
'ball'
: the spatial parameter is the radius , with
, the temporal parameter is the temporal lag
:
'gauss'
: Gaussian function where 95% of the mass is inside the spatial parameter , with
, the temporal parameter is the temporal lag
:
Above, represents the indicator function. The argument
kernel_type
determines the used kernel function as presented above, the argument lags
provides the used temporal lags for the kernel functions ( in the above formulas) and the argument
kernel_parameters
gives the spatial parameters for the kernel function. Each of the arguments kernel_type
, lags
and kernel_parameters
can be of length K
or 1
. Specifically, kernel_type
can be either one kernel, then each local autocovariance matrix use the same kernel type, or of length K
which leads to different kernel functions for the provided kernel parameters. lags
can be either one integer, then for each kernel the same temporal lag is used, or an integer vector of length K
which leads to different temporal lags. In the same fashion kernel_parameters
is a vector of length K
or 1
. If kernel_type
equals 'ball'
or 'gauss'
then the corresponding entry of kernel_parameters
gives the single spatial radius parameter. In contrast, if (at least one entry of) kernel_type
equals 'ring'
, then kernel_parameters
must be a list of length K
(or 1
) where each entry is a numeric vector of length 2
defining the inner and outer spatial radius.
Internally, stbss
calls stkmat
to compute a list of c(n,n)
kernel matrices based on the parameters given, where each entry of those matrices corresponds to . Alternatively, such a list of kernel matrices can be given directly to the function
stbss
through the kernel_list
argument. This is useful when stbss
is called numerous times with the same coordinates/kernel functions as the computation of the kernel matrices is then only done once prior the actual stbss
calls. For details see also lacov
.
If more than one local autocovariance matrix is used stbss
jointly diagonalizes these matrices with the function frjd
. ...
provides arguments for frjd
, useful arguments might be:
eps
: tolerance for convergence.
maxiter
: maximum number of iterations.
stbss
returns a list of class 'stbss'
with the following entries:
s |
object of |
coords |
coordinates of the observations. Is |
time |
time of the observations. Is |
w |
estimated unmixing matrix. |
w_inv |
inverse of the estimated unmixing matrix. |
pevals |
(pseudo-)eigenvalues for each latent field entry. |
d |
matrix of stacked (jointly) diagonalized local autocovariance matrices with dimension |
x_mu |
columnmeans of |
cov_inv_sqrt |
square root of the inverse sample covariance matrix of |
Muehlmann, C., De Iaco, S. and Nordhausen, K. (2023), Blind Recovery of Sources for Multivariate Space-Time Environmental Data. Stochastic and Environmental Research and Risk Assessment, 37, 1593–1613, <doi:10.1007/s00477-022-02348-2>.
# space and time coordinates n_t <- 50 n_sp <- 10 st_coords <- as.matrix(expand.grid(1:n_sp, 1:n_sp, 1:n_t)) # simulate three latent white noise fields field_1 <- rnorm(nrow(st_coords)) field_2 <- rnorm(nrow(st_coords)) field_3 <- rnorm(nrow(st_coords)) # compute the observed field latent_field <- cbind(field_1, field_2, field_3) mixing_matrix <- matrix(rnorm(9), 3, 3) observed_field <- latent_field # apply stbss with lag 1 and a ring kernel stbss_res <- stbss(observed_field, coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ring', kernel_parameters = list(c(0, 1)), lags = 1) # print object print(stbss_res) # unmixing matrix w_unmix <- coef(stbss_res) # apply the same stbss with a kernel list kernel_list <- stkmat(coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ring', kernel_parameters = list(c(0, 1)), lags = 1) stbss_res_k <- stbss(observed_field, kernel_list = kernel_list) # apply stbss with three ball kernels stbss_res_b <- stbss(observed_field, coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ball', kernel_parameters = 1:3, lags = 1:3)
# space and time coordinates n_t <- 50 n_sp <- 10 st_coords <- as.matrix(expand.grid(1:n_sp, 1:n_sp, 1:n_t)) # simulate three latent white noise fields field_1 <- rnorm(nrow(st_coords)) field_2 <- rnorm(nrow(st_coords)) field_3 <- rnorm(nrow(st_coords)) # compute the observed field latent_field <- cbind(field_1, field_2, field_3) mixing_matrix <- matrix(rnorm(9), 3, 3) observed_field <- latent_field # apply stbss with lag 1 and a ring kernel stbss_res <- stbss(observed_field, coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ring', kernel_parameters = list(c(0, 1)), lags = 1) # print object print(stbss_res) # unmixing matrix w_unmix <- coef(stbss_res) # apply the same stbss with a kernel list kernel_list <- stkmat(coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ring', kernel_parameters = list(c(0, 1)), lags = 1) stbss_res_k <- stbss(observed_field, kernel_list = kernel_list) # apply stbss with three ball kernels stbss_res_b <- stbss(observed_field, coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ball', kernel_parameters = 1:3, lags = 1:3)
Computation of spatio-temporal kernel matrices for given kernel functions.
stkmat(coords, time, kernel_type, kernel_parameters, lags)
stkmat(coords, time, kernel_type, kernel_parameters, lags)
coords |
a numeric matrix of dimension |
time |
an integer vector of length |
kernel_type |
either a string or a string vector of length |
kernel_parameters |
a numeric vector of length |
lags |
an integer vector of length |
The following kernel functions are implemented and chosen with the argument kernel_type
:
'ring'
: the spatial parameters are inner radius and outer radius
, with
, and
, the temporal parameter is the temporal lag
:
'ball'
: the spatial parameter is the radius , with
, the temporal parameter is the temporal lag
:
'gauss'
: Gaussian function where 95% of the mass is inside the spatial parameter , with
, the temporal parameter is the temporal lag
:
Above, represents the indicator function. The argument
kernel_type
determines the used kernel function as presented above, the argument lags
provides the used temporal lags for the kernel functions ( in the above formulas) and the argument
kernel_parameters
gives the spatial parameters for the kernel function. Each of the arguments kernel_type
, lags
and kernel_parameters
can be of length K
or 1
. Specifically, kernel_type
can be either one kernel, then each local autocovariance matrix use the same kernel type, or of length K
which leads to different kernel functions for the provided kernel parameters. lags
can be either one integer, then for each kernel the same temporal lag is used, or an integer vector of length K
which leads to different temporal lags. In the same fashion kernel_parameters
is a vector of length K
or 1
. If kernel_type
equals 'ball'
or 'gauss'
then the corresponding entry of kernel_parameters
gives the single spatial radius parameter. In contrast, if (at least one entry of) kernel_type
equals 'ring'
, then kernel_parameters
must be a list of length K
(or 1
) where each entry is a numeric vector of length 2
defining the inner and outer spatial radius. See examples below.
The output of this function can be used with the function stbss
to avoid unnecessary computation of kernel matrices when stbss
is called multiple times with the same coordinate/kernel function setting. Additionally, the output can be used with the function lacov
.
stkmat
returns a list of length K
containing numeric matrices of dimension c(n,n)
corresponding to the spatio-temporal kernel matrices.
Muehlmann, C., De Iaco, S. and Nordhausen, K. (2023), Blind Recovery of Sources for Multivariate Space-Time Environmental Data. Stochastic and Environmental Research and Risk Assessment, 37, 1593–1613, <doi:10.1007/s00477-022-02348-2>.
# space and time coordinates n_t <- 50 n_sp <- 10 coords <- runif(n_sp ^ 2 * 2) * n_sp dim(coords) <- c(n_sp ^ 2, 2) time <- 1:n_t st_coords <- as.matrix(expand.grid(1:nrow(coords), 1:length(time))) st_coords <- cbind(coords[st_coords[, 1], ], time[st_coords[, 2]]) # different ring kernels and same lags stkmat_r <- stkmat(coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ring', kernel_parameters = list(c(0, 1), c(1, 2)), lags = c(1, 1)) # same ball kernels and different lags stkmat_b <- stkmat(coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ball', kernel_parameters = 1:3, lags = c(1, 2, 3)) # different gauss kernels and different lags stkmat_g <- stkmat(coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'gauss', kernel_parameters = 1:3, lags = 1:3) # mixed kernels stkmat_m <- stkmat(coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = c('ball', 'ring', 'gauss'), kernel_parameters = list(1, c(1:2), 3), lags = 1:3)
# space and time coordinates n_t <- 50 n_sp <- 10 coords <- runif(n_sp ^ 2 * 2) * n_sp dim(coords) <- c(n_sp ^ 2, 2) time <- 1:n_t st_coords <- as.matrix(expand.grid(1:nrow(coords), 1:length(time))) st_coords <- cbind(coords[st_coords[, 1], ], time[st_coords[, 2]]) # different ring kernels and same lags stkmat_r <- stkmat(coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ring', kernel_parameters = list(c(0, 1), c(1, 2)), lags = c(1, 1)) # same ball kernels and different lags stkmat_b <- stkmat(coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'ball', kernel_parameters = 1:3, lags = c(1, 2, 3)) # different gauss kernels and different lags stkmat_g <- stkmat(coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = 'gauss', kernel_parameters = 1:3, lags = 1:3) # mixed kernels stkmat_m <- stkmat(coords = st_coords[, 1:2], time = st_coords[, 3], kernel_type = c('ball', 'ring', 'gauss'), kernel_parameters = list(1, c(1:2), 3), lags = 1:3)