Merge branch 'stable' into cran

.Rbuildignore

@@ -13,9 +13,4 @@
 ^shiny_app$
 ^.vscode$
 ^deprecated$
-^vignettes/.*\.Rmd$
-^vignettes/.*\.Rmd$
-^vignettes/.*\.Rmd$
-^vignettes/.*\.Rmd$
-^vignettes/.*\.Rmd$
-^vignettes/.*\.Rmd$
+^cran_vignettes$

.gitignore

@@ -34,6 +34,7 @@ docs
 
 # produced vignettes
 vignettes/*.pdf
+cran_vignettes/*.html
 
 # OAuth2 token, see https://github.com/hadley/httr/releases/tag/v0.3
 .httr-oauth

DESCRIPTION

@@ -2,22 +2,22 @@ Package: rSPDE
 Type: Package
 Title: Rational Approximations of Fractional Stochastic Partial
         Differential Equations
-Version: 2.3.2
+Version: 2.3.3
 Authors@R: c(
     person("David", "Bolin", email = "[email protected]", role = c("cre", "aut")),
     person("Alexandre", "Simas", email = "[email protected]", role = "aut"),
     person("Finn", "Lindgren", email = "[email protected]", role = "ctb"))
 Maintainer: David Bolin <[email protected]>
 Description: Functions that compute rational approximations of fractional elliptic stochastic partial differential equations. The package also contains functions for common statistical usage of these approximations. The main references for rSPDE are Bolin, Simas and Xiong (2023) <doi:10.1080/10618600.2023.2231051> for the covariance-based method and Bolin and Kirchner (2020) <doi:10.1080/10618600.2019.1665537> for the operator-based rational approximation. These can be generated by the citation function in R.
 Depends: R (>= 3.5.0), Matrix
-Imports: stats, methods
+Imports: stats, methods, fmesher, lifecycle, broom
 License: GPL (>=3) | file LICENSE
 URL: https://davidbolin.github.io/rSPDE/
 Encoding: UTF-8
 RoxygenNote: 7.2.3
-Suggests: R.rsp, rmarkdown, INLA (>= 22.12.14), testthat, rgdal,
+Suggests: knitr, rmarkdown, INLA (>= 22.12.14), testthat,
         ggplot2, lattice, splancs, optimParallel, RSpectra, numDeriv,
-        inlabru (>= 2.7.0), sn, viridis, scoringRules, doParallel, foreach
+        inlabru (>= 2.7.0), sn, viridis, scoringRules, doParallel, foreach, tidyr
 Additional_repositories: https://inla.r-inla-download.org/R/testing
 BugReports: https://github.com/davidbolin/rSPDE/issues
 VignetteBuilder: R.rsp

NAMESPACE

@@ -1,7 +1,11 @@
 # Generated by roxygen2: do not edit by hand
 
+S3method(augment,rspde_lme)
+S3method(deviance,rspde_lme)
 S3method(gg_df,rspde_result)
+S3method(glance,rspde_lme)
 S3method(logLik,rspde_lme)
+S3method(nobs,rspde_lme)
 S3method(precision,CBrSPDEobj)
 S3method(precision,inla_rspde)
 S3method(predict,CBrSPDEobj)
@@ -36,20 +40,24 @@ export(Qsqrt.mult)
 export(Qsqrt.solve)
 export(Sigma.mult)
 export(Sigma.solve)
+export(augment)
 export(construct.spde.matern.loglike)
 export(cross_validation)
 export(folded.matern.covariance.1d)
 export(folded.matern.covariance.2d)
 export(fractional.operators)
 export(get.initial.values.rSPDE)
 export(gg_df)
+export(glance)
+export(graph_data_rspde)
 export(intrinsic.matern.operators)
 export(matern.covariance)
 export(matern.operators)
 export(precision)
 export(rSPDE.A1d)
 export(rSPDE.construct.matern.loglike)
 export(rSPDE.fem1d)
+export(rSPDE.fem2d)
 export(rSPDE.loglike)
 export(rSPDE.matern.loglike)
 export(rational.order)
@@ -67,6 +75,7 @@ export(rspde.mesh.projector)
 export(rspde.metric_graph)
 export(rspde.result)
 export(rspde_lme)
+export(spde.make.A)
 export(spde.matern.loglike)
 export(spde.matern.operators)
 export(variogram.intrinsic.spde)
@@ -78,13 +87,24 @@ if (getRversion() >= "3.6.0") {
 }
 
 import(Matrix)
+importFrom(broom,augment)
+importFrom(broom,glance)
+importFrom(fmesher,fm_basis)
+importFrom(fmesher,fm_block)
+importFrom(fmesher,fm_fem)
+importFrom(fmesher,fm_mesh_1d)
+importFrom(fmesher,fm_mesh_2d)
+importFrom(fmesher,fm_row_kron)
+importFrom(lifecycle,deprecated)
 importFrom(methods,as)
 importFrom(stats,approx)
 importFrom(stats,dbeta)
+importFrom(stats,deviance)
 importFrom(stats,dnorm)
 importFrom(stats,lm)
 importFrom(stats,logLik)
 importFrom(stats,na.omit)
+importFrom(stats,nobs)
 importFrom(stats,optim)
 importFrom(stats,pnorm)
 importFrom(stats,quantile)

NEWS.md

@@ -1,3 +1,14 @@
+# rSPDE 2.3.3
+* Bugfix on rspde_lme when fitting with fixed smoothness.
+* Added a 2d fem interface.
+* Moved from using INLA's mesh functions to fmesher's mesh functions.
+* Removing rgdal from suggests.
+* The `data` argument in `predict.rspde_lme` has been changed to `newdata`.
+* Adding `covariance_mesh` and `cov_function_mesh` methods as functions in the list returned by objects obtained from `matern.operators()` and `spde.matern.operators()`.
+* Updated the internal structure to match the updates from the `MetricGraph` package.
+* Updated the `cross_validation` function to match the updates in `inlabru`.
+* Added `glance` and `augment` methods for `rspde_lme` objects.
+
 # rSPDE 2.3.2
 * Small improvement on speed for rspde_lme.
 * Bugfix on Q for small values of nu in dimension 1.

R/fractional.computations.R

@@ -259,28 +259,48 @@ update.CBrSPDEobj <- function(object, user_nu = NULL, user_alpha = NULL,
           graph <- new_object$graph
         }
 
-        new_object <- matern.operators(
-          kappa = new_object$kappa,
-          sigma = new_object$sigma,
-          range = new_object$range,
-          tau = new_object$tau,
-          nu = new_object$nu,
-          alpha = new_object$alpha,
-          G = new_object$G,
-          C = new_object$C,
-          d = new_object$d,
-          m = new_object$m,
-          mesh = mesh,
-          loc_mesh = loc_mesh,
-          range_mesh = range_mesh,
-          graph = graph,
-          parameterization = parameterization,
-          type = "covariance",
-          return_block_list = return_block_list,
-          type_rational_approximation = type_rational_approximation,
-          fem_mesh_matrices = new_object$fem_mesh_matrices,
-          compute_logdet = new_object$compute_logdet
-        )
+        if(parameterization == "spde"){
+            new_object <- matern.operators(
+              kappa = new_object$kappa,
+              tau = new_object$tau,
+              alpha = new_object$alpha,
+              G = new_object$G,
+              C = new_object$C,
+              d = new_object$d,
+              m = new_object$m,
+              mesh = mesh,
+              loc_mesh = loc_mesh,
+              range_mesh = range_mesh,
+              graph = graph,
+              parameterization = parameterization,
+              type = "covariance",
+              return_block_list = return_block_list,
+              type_rational_approximation = type_rational_approximation,
+              fem_mesh_matrices = new_object$fem_mesh_matrices,
+              compute_logdet = new_object$compute_logdet
+            )
+        } else{
+            new_object <- matern.operators(
+              sigma = new_object$sigma,
+              range = new_object$range,
+              nu = new_object$nu,
+              G = new_object$G,
+              C = new_object$C,
+              d = new_object$d,
+              m = new_object$m,
+              mesh = mesh,
+              loc_mesh = loc_mesh,
+              range_mesh = range_mesh,
+              graph = graph,
+              parameterization = parameterization,
+              type = "covariance",
+              return_block_list = return_block_list,
+              type_rational_approximation = type_rational_approximation,
+              fem_mesh_matrices = new_object$fem_mesh_matrices,
+              compute_logdet = new_object$compute_logdet
+            )
+        }
+
   } else{
   ## get parameters
 
@@ -340,12 +360,11 @@ update.CBrSPDEobj <- function(object, user_nu = NULL, user_alpha = NULL,
           alpha <- new_object$nu + d / 2
           new_object$alpha <- alpha
         } 
-
+        if(parameterization == "spde"){
           new_object <- spde.matern.operators(
             kappa = new_object$kappa,
             tau = new_object$tau,
             theta = new_object$theta,
-            nu = new_object$nu,
             alpha = new_object$alpha,
             B.tau = new_object$B.tau,
             B.kappa = new_object$B.kappa,
@@ -358,11 +377,31 @@ update.CBrSPDEobj <- function(object, user_nu = NULL, user_alpha = NULL,
             range_mesh = range_mesh,
             graph = graph,            
             parameterization = parameterization,
+            type = "covariance",
+            type_rational_approximation = new_object$type_rational_approximation
+          )
+        } else{
+          new_object <- spde.matern.operators(
+            kappa = new_object$kappa,
+            tau = new_object$tau,
+            theta = new_object$theta,
+            nu = new_object$nu,
+            G = new_object$G,
+            C = new_object$C,
+            d = new_object$d,
+            m = new_object$m,
+            mesh = mesh,
+            loc_mesh = loc_mesh,
+            range_mesh = range_mesh,
+            graph = graph,            
+            parameterization = parameterization,
             B.sigma = new_object$B.sigma,
             B.range = new_object$B.range,
             type = "covariance",
             type_rational_approximation = new_object$type_rational_approximation
           )
+        }
+
   }
   return(new_object)
 }
@@ -513,24 +552,40 @@ update.rSPDEobj <- function(object, user_nu = NULL,
           graph <- new_object$graph
         }
 
-  new_object <- matern.operators(
-    kappa = new_object$kappa,
-    sigma = new_object$sigma,
-    range = new_object$range,
-    tau = new_object$tau,
-    nu = new_object$nu,
-    alpha = new_object$alpha,
-    G = new_object$G,
-    C = new_object$C,
-    d = new_object$d,
-    m = new_object$m,
-    mesh = mesh,
-    loc_mesh = loc_mesh,
-    range_mesh = range_mesh,
-    graph = graph,
-    parameterization = parameterization,    
-    type = "operator"
-  )
+        if(parameterization == "spde"){
+            new_object <- matern.operators(
+              kappa = new_object$kappa,
+              tau = new_object$tau,
+              alpha = new_object$alpha,
+              G = new_object$G,
+              C = new_object$C,
+              d = new_object$d,
+              m = new_object$m,
+              mesh = mesh,
+              loc_mesh = loc_mesh,
+              range_mesh = range_mesh,
+              graph = graph,
+              parameterization = parameterization,    
+              type = "operator"
+            )
+        } else{
+            new_object <- matern.operators(
+              sigma = new_object$sigma,
+              range = new_object$range,
+              nu = new_object$nu,
+              G = new_object$G,
+              C = new_object$C,
+              d = new_object$d,
+              m = new_object$m,
+              mesh = mesh,
+              loc_mesh = loc_mesh,
+              range_mesh = range_mesh,
+              graph = graph,
+              parameterization = parameterization,    
+              type = "operator"
+            )
+        }
+
  } else{
           if (!is.null(user_tau)) {
             new_object$tau <- rspde_check_user_input(user_tau, "tau", 0)
@@ -590,12 +645,12 @@ update.rSPDEobj <- function(object, user_nu = NULL,
 
         }         
 
-
+    if(parameterization == "spde"){
           new_object <- spde.matern.operators(
             kappa = new_object$kappa,
             tau = new_object$tau,
             theta = new_object$theta,
-            nu = new_object$nu,
+            alpha = new_object$alpha,
             B.tau = new_object$B.tau,
             B.kappa = new_object$B.kappa,
             G = new_object$G,
@@ -609,6 +664,26 @@ update.rSPDEobj <- function(object, user_nu = NULL,
             graph = graph,            
             type = "operator"
           )  
+    } else{
+        new_object <- spde.matern.operators(
+            kappa = new_object$kappa,
+            tau = new_object$tau,
+            theta = new_object$theta,
+            nu = new_object$nu,
+            B.range = new_object$B.range,
+            B.sigma = new_object$B.sigma,
+            G = new_object$G,
+            C = new_object$C,
+            d = new_object$d,
+            m = new_object$m,
+            mesh = mesh,
+            loc_mesh = loc_mesh,
+            range_mesh = range_mesh,
+            parameterization = parameterization,
+            graph = graph,            
+            type = "operator"
+          )  
+    }
  }
 
 
@@ -712,14 +787,18 @@ simulate.CBrSPDEobj <- function(object, nsim = 1,
     sizeL <- dim(L)[1]
 
 
-    Q <- rspde.matern.precision.integer(
-      kappa = kappa, nu = nu, tau = tau,
-      dim = d,
-      fem_mesh_matrices = fem_mesh_matrices
-    )
+    # Q <- rspde.matern.precision.integer(
+    #   kappa = kappa, nu = nu, tau = tau,
+    #   dim = d,
+    #   fem_mesh_matrices = fem_mesh_matrices
+    # )
+
+    Q <- object$Q
+
     Z <- rnorm(sizeL * nsim)
     dim(Z) <- c(sizeL, nsim)
 
+
     # LQ <- Matrix::Cholesky(forceSymmetric(Q), LDL = FALSE)
     # X <- solve(LQ, Z, system = "Lt")
     # X <- solve(LQ, X, system = "Pt")
@@ -963,7 +1042,7 @@ predict.rSPDEobj <- function(object,
 #' op <- matern.operators(
 #'   range = range, sigma = sigma, nu = nu,
 #'   loc_mesh = x, d = 1,
-#'   type = "operator"
+#'   type = "operator", parameterization = "matern"
 #' )
 #'
 #' # Sample the model