View source: R/smoothing_time.R
smooth.FEM.time | R Documentation |
Space-time regression with differential regularization. Space-varying covariates can be included in the model. The technique accurately handle data distributed over irregularly shaped domains. Moreover, various conditions can be imposed at the domain boundaries.
smooth.FEM.time(locations = NULL, time_locations = NULL, observations, FEMbasis,
time_mesh=NULL, covariates = NULL, PDE_parameters = NULL, BC = NULL,
incidence_matrix = NULL, areal.data.avg = TRUE,
FLAG_MASS = FALSE, FLAG_PARABOLIC = FALSE, FLAG_ITERATIVE = FALSE,
threshold = 10^(-4), max.steps = 50, IC = NULL,
search = "tree", bary.locations = NULL,
family = "gaussian", mu0 = NULL, scale.param = NULL,
threshold.FPIRLS = 0.0002020, max.steps.FPIRLS = 15,
lambda.selection.criterion = "grid", DOF.evaluation = NULL,
lambda.selection.lossfunction = NULL, lambdaS = NULL, lambdaT = NULL,
DOF.stochastic.realizations = 100, DOF.stochastic.seed = 0,
DOF.matrix = NULL, GCV.inflation.factor = 1, lambda.optimization.tolerance = 0.05,
inference.data.object.time=NULL)
locations |
A matrix where each row specifies the spatial coordinates |
time_locations |
A vector containing the times of the corresponding observations in the vector |
observations |
A matrix of #locations x #time_locations with the observed data values over the spatio-temporal domain.
The spatial locations of the observations can be specified with the |
FEMbasis |
A |
time_mesh |
A vector specifying the time mesh. |
covariates |
A #observations-by-#covariates matrix where each row represents the covariates associated with the corresponding observed data value in |
PDE_parameters |
A list specifying the parameters of the PDE in the regularizing term. Default is NULL, i.e. regularization is by means of the Laplacian (stationary, isotropic case).
If the PDE is elliptic it must contain: |
BC |
A list with two vectors:
|
incidence_matrix |
A #regions-by-#triangles/tetrahedrons matrix where the element (i,j) equals 1 if the j-th triangle/tetrahedron is in the i-th region and 0 otherwise.
This is only for areal data. In case of pointwise data, this parameter is set to |
areal.data.avg |
Boolean. It involves the computation of Areal Data. If |
FLAG_MASS |
Boolean. This parameter is considered only for separable problems i.e. when |
FLAG_PARABOLIC |
Boolean. If |
FLAG_ITERATIVE |
Boolean. If |
threshold |
This parameter is used for arresting algorithm iterations. Algorithm stops when two successive iterations lead to improvement in penalized log-likelihood smaller than threshold.
Default value |
max.steps |
This parameter is used to limit the maximum number of iteration.
Default value |
IC |
Initial condition needed in case of parabolic problem i.e. when |
search |
a flag to decide the search algorithm type (tree or naive or walking search algorithm). |
bary.locations |
A list with three vectors:
|
family |
This parameter specify the distribution within exponential family used for GLM model.
The following distribution are implemented: "binomial", "exponential", "gamma", "poisson", "gaussian", "invgaussian".
The default link function for binomial is |
mu0 |
This parameter is a vector that set the starting point for FPIRLS algorithm. It represent an initial guess of the location parameter.
Default is set to observation for non binary distribution while equal to |
scale.param |
Dispersion parameter of the chosen distribution. This is only required for "gamma", "gaussian", "invgaussian". User may specify the parameter as a positive real number. If the parameter is not supplied, it is estimated from data according to Wilhelm Sangalli 2016. |
threshold.FPIRLS |
This parameter is used for arresting algorithm iterations. Algorithm stops when two successive iterations lead to improvement in penalized log-likelihood smaller than threshold.FPIRLS.
Default value |
max.steps.FPIRLS |
This parameter is used to limit the maximum number of iteration.
Default value |
lambda.selection.criterion |
This parameter is used to select the optimization method related to smoothing parameter |
DOF.evaluation |
This parameter is used to identify if and how degrees of freedom computation has to be performed.
The following possibilities are allowed: NULL, 'exact' and 'stochastic'
In the former case no degree of freedom is computed, while the other two methods enable computation.
Stochastic computation of DOFs may be slightly less accurate than its deterministic counterpart, but is highly suggested for meshes of more than 5000 nodes, being fairly less time consuming.
Default value |
lambda.selection.lossfunction |
This parameter is used to understand if some loss function has to be evaluated.
The following possibilities are allowed: NULL and 'GCV' (generalized cross validation)
The former case is that of |
lambdaS |
A scalar or vector of spatial smoothing parameters. |
lambdaT |
A scalar or vector of temporal smoothing parameters. |
DOF.stochastic.realizations |
This parameter is considered only when |
DOF.stochastic.seed |
This parameter is considered only when |
DOF.matrix |
Matrix of degrees of freedom. This parameter can be used if the DOF.matrix corresponding to |
GCV.inflation.factor |
Tuning parameter used for the estimation of GCV. Default value |
lambda.optimization.tolerance |
Tolerance parameter, a double between 0 and 1 that fixes how much precision is required by the optimization method: the smaller the parameter, the higher the accuracy.
Used only if |
inference.data.object.time |
An |
A list with the following variables:
fit.FEM.time
A FEM.time
object that represents the fitted spatio-temporal field.
PDEmisfit.FEM.time
A FEM.time
object that represents the misfit of the penalized PDE.
beta
If covariates
is not NULL
, a matrix with number of rows equal to the number of covariates and number of columns equal to length of lambda. The j
th column represents the vector of regression coefficients when
the smoothing parameter is equal to lambda[j]
.
edf
If GCV is TRUE
, a scalar or matrix with the trace of the smoothing matrix for each combination of the smoothing parameters specified in lambdaS
and lambdaT
.
stderr
If GCV is TRUE
, a scalar or matrix with the estimate of the standard deviation of the error for each combination of the smoothing parameters specified in lambdaS
and lambdaT
.
GCV
If GCV is TRUE
, a scalar or matrix with the value of the GCV criterion for each combination of the smoothing parameters specified in lambdaS
and lambdaT
.
bestlambda
If GCV is TRUE
, a 2-elements vector with the indices of smoothing parameters returning the lowest GCV
ICestimated
If FLAG_PARABOLIC is TRUE
and IC is NULL
, a list containing a FEM
object with the initial conditions, the value of the smoothing parameter lambda returning the lowest GCV and, in presence of covariates, the estimated beta coefficients
bary.locations
A barycenter information of the given locations if the locations are not mesh nodes.
inference
A list set only if a well defined inferenceDataObjectTime
is passed as parameter to the function; contains all inference outputs required:
p_values
list of lists set only if at least one p-value is required; contains the p-values divided by implementation:
wald
list containing all the Wald p-values required, in the same order of the type
list in inference.data.object.time
.
If one-at-the-time tests are required, the corresponding item is a vector of p values ordered as the rows of coeff
matrix in inference.data.object.time
.
speckman
list containing all the Speckman p-values required, in the same order of the type
list in inference.data.object.time
.
If one-at-the-time tests are required, the corresponding item is a vector of p values ordered as the rows of coeff
matrix in inference.data.object.time
.
eigen_sign_flip
list containing all the Eigen-Sign-Flip p-values required, in the same order of the type
list in inference.data.object.time
.
If one-at-the-time tests are required, the corresponding item is a vector of p values ordered as the rows of coeff
matrix in inference.data.object.time
.
CI
list of lists set only if at least one confidence interval is required; contains the confidence intervals divided by implementation:
wald
list containing all the Wald confidence intervals required, in the same order of the type
list in inference.data.object.time
.
Each item is a matrix with 3 columns and p rows, p being the number of rows of coeff
matrix in inference.data.object.time
; each row is the CI for the corresponding row of coeff
matrix.
speckman
list containing all the Speckman confidence intervals required, in the same order of the type
list in inference.data.object.time
.
Each item is a matrix with 3 columns and p rows, p being the number of rows of coeff
matrix in inference.data.object.time
; each row is the CI for the corresponding row of coeff
matrix.
#' @references Arnone, E., Azzimonti, L., Nobile, F., & Sangalli, L. M. (2019). Modeling spatially dependent functional data via regression with differential regularization. Journal of Multivariate Analysis, 170, 275-295. Bernardi, M. S., Sangalli, L. M., Mazza, G., & Ramsay, J. O. (2017). A penalized regression model for spatial functional data with application to the analysis of the production of waste in Venice province. Stochastic Environmental Research and Risk Assessment, 31(1), 23-38. Federico Ferraccioli, Laura M. Sangalli & Livio Finos (2022). Some first inferential tools for spatial regression with differential regularization. Journal of Multivariate Analysis, 189, 104866.
library(fdaPDE)
data(horseshoe2D)
boundary_nodes = horseshoe2D$boundary_nodes
boundary_segments = horseshoe2D$boundary_segments
locations = horseshoe2D$locations
time_locations = seq(0,1,length.out = 5)
mesh = create.mesh.2D(nodes = rbind(boundary_nodes, locations), segments = boundary_segments)
space_time_locations = cbind(rep(time_locations,each=nrow(mesh$nodes)),
rep(mesh$nodes[,1],5),rep(mesh$nodes[,2],5))
FEMbasis = create.FEM.basis(mesh)
lambdaS = 10^-1
lambdaT = 10^-1
data = fs.test(space_time_locations[,2],
space_time_locations[,3])*cos(pi*space_time_locations[,1]) +
rnorm(nrow(space_time_locations), sd = 0.5)
data = matrix(data, nrow = nrow(mesh$nodes), ncol = length(time_locations), byrow = TRUE)
solution = smooth.FEM.time(observations = data, time_locations = time_locations,
FEMbasis = FEMbasis, lambdaS = lambdaS, lambdaT = lambdaT)
plot(solution$fit.FEM)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.