R/smooth.construct.dm.smooth.spec.R

In dill/msg: Multidimensional Scaling for Generalized Additive Models

#' Distance matrix smoothing
#'
#' Provides smoothing methods for multidimensional scaling-based projections of a given distance matrix. The matrix can be supplied or a function to calculate the distances can be supplied.
#'
#' Smoothing is performed using Duchon splines (see \code{\link{Duchon.spline}} for more information).
#'
#' @aliases dm smooth.construct.dm.smooth.spec
#' @method smooth.construct dm.smooth.spec
#' @export
#' @import mgcv
#' @useDynLib msg wood_path
#'
#' @param object a smooth specification object, usually generated by a term \code{s(...,bs="ds",...)}. Note that \code{xt} object is needed, see Details.
#' @param data a list containing just the data (including any \code{by} variable) required by this term, with names corresponding to \code{object$term} (and \code{object$by}). The \code{by} variable is the last element.
#' @param knots IGNORED!
#'
#' @return An object of class \code{dm.smooth}. In addition to the usual elements of a smooth class documented under \code{\link{smooth.construct}}, this object will contain an element named \code{msg}:
#' \describe{
#'   \item{\code{mds.obj}}{result of running \code{\link{cmdscale}} on the data.}
#'   \item{\code{dim}}{dimension of the MDS projection.}
#'   \item{\code{term}}{auto-generated names of the variables in the MDS space (of the form "mds-i" where i indexes the data).}
#'   \item{\code{data}}{the data projected into MDS space.}
#'   \item{\code{\dots}}{Plus those extra elements as documented in \code{\link{Duchon.spline}}}}
#'
#' @section Details: The constructor is not normally called directly, but is rather used internally by \code{\link{gam}}. To use for basis setup it is recommended to use \code{\link{smooth.construct2}}.
#'
#' When specifying the model extra arguments must be supplied by the \code{xt} argument.
#' \describe{
#'   \item{\code{D}}{a distance matrix}
#'   \item{\code{mds.dim}}{dimension of the MDS projection}
#'   \item{\code{grid}}{a grid over the covariates to use as a base for the MDS configuration. If \code{NULL} the sample points will be used. See below for details.}
#'   \item{\code{dist_fn}}{function to calculate distances, see below}}
#'
#' **Write something here about grid**
#'
#' \code{dist_fn} takes one argument, a \code{data.frame} of locations of the data, as provided as the covariates used in the smooth.
#'
#' MDS dimension selection may be performed by finding the projection with the lowest GCV score. BEWARE: the GCV score is not necessarily monotonic in the number of dimensions. Automated dimension selection will appear in a later version of the package.
#'
#' @references
#' Duchon, J. (1977) Splines minimizing rotation-invariant semi-norms in Solobev spaces. in W. Shemp and K. Zeller (eds) Construction theory of functions of several variables, 85-100, Springer, Berlin.
#' Miller, DL and Wood, SN. (2014) Finite area smoothing with generalized distance splines. Environmental and Ecological Statistics 4 715-731
#'
#' @author David L. Miller
#'
#' @examples
#' \dontrun{
#' # test this works with the wt2 example from msg
#' library(msg)
#' data(wt2)
#'
#' ## using a pre-built D matrix
#' # create the sample
#' samp.ind <- sample(1:length(wt2$data$x),250)
#' wt2.samp <- list(x=wt2$data$x[samp.ind],
#'                 y=wt2$data$y[samp.ind],
#'                 z=wt2$data$z[samp.ind]+rnorm(250)*0.9)
#' mds.dim<-5
#' custom_dist_fn <- function(x){
#'   msg:::create_distance_matrix(x$x,x$y,wt2$bnd,faster=0)
#' }
#'
#' grid_obj <- msg:::create_refgrid(wt2$bnd,120)
#'
#' grid_obj$nrefx <- grid_obj$nrefy <- NULL
#' grid_obj <- as.data.frame(grid_obj)
#'
#'
#' b.dm <- gam(z~s(x, y, bs="dm", k=200,
#'                 xt=list(dist_fn=custom_dist_fn, mds.dim=5, grid=grid_obj)),
#'             data=wt2.samp)
#'
#'
#' # with msg
#' b.msg<-gam(z~s(x,y,bs="msg",k=200,xt=list(bnd=wt2$bnd,mds.dim=5)),data=wt2.samp)
#' }
#' @keywords models smoothing
smooth.construct.dm.smooth.spec <- function(object, data, knots){

   # this is the dm smooth.spec
   # this does most of the work
   if(is.null(object$xt$mds.dim)){
      stop("No MDS projection dimension supplied!\n")
   }

   # extract the MDS dimension
   mds.dim <- object$xt$mds.dim

   new.obj <- list()

   # use the distance function
   if(is.null(object$xt$D) & !is.null(object$xt$dist_fn)){
     dist_fn <- object$xt$dist_fn
     if(!is.null(object$xt$grid)){
       grid_points <- object$xt$grid
       # get the grid distances
       D_grid <- dist_fn(grid_points)
       # MDS them
       mds_obj <- cmdscale(D_grid, eig=TRUE, k=mds.dim, x.ret=TRUE)
       # insert the sample in
       mds_sample <- insert.mds.generic(mds_obj, data, grid_points, dist_fn)

     # get grid distances
     }else{
     # if there is no grid then just use the samples
       D_sample <- dist_fn(data)
       # MDS them
       mds_obj <- cmdscale(D_grid, eig=TRUE, k=mds.dim, x.ret=TRUE)

       mds_sample <- mds_obj$points
    }

   }else if(!is.null(object$xt$D) & is.null(object$xt$dist_fn)){
   # use the saved distance matrix
      new.obj$D <- object$xt$D
      D <- object$xt$D

      # do the mds projection
      mds_obj <- cmdscale(D, eig=TRUE, k=mds.dim, x.ret=TRUE)
      mds_sample <- mds.obj$points
   }else{
     stop("Need to specify either a distance matrix or distance function")
   }

   # save the mds stuff
   object$xt$mds <- mds_obj

   # make some variable names up
   mds.names <- paste("mds-",1:dim(mds_sample)[2],sep="")
   # make sure the data is a data.frame
   mds_sample <- as.data.frame(mds_sample)
   # remove any already in the data
   names(mds_sample) <- mds.names

   # make sure there are the right stuff is in the object before passing
   # to Duchon, but save beforehand!
   save.dim <- object$dim
   save.term <- object$term

   # reset the object elements
   object$term <- mds.names
   object$dim <- mds.dim
   #object$xt$mds_data <- mds_sample

   # if knots were supplied, they're going to be ignored, warn about that!
   if(length(knots)!=0){
      warning("Knots were supplied but will be ignored!\n")
      knots <- list()
   }

   # set the penalty order
   object$p.order <- c(2,mds.dim/2-1)

   # make the duchon splines object as usual
   object <- smooth.construct.ds.smooth.spec(object, mds_sample, knots)

   # recover the stuff we want in the object
   object$term <- save.term
   object$dim <- save.dim

   class(object) <- "dm.smooth"
   object
}


#' Distance matrix smoothing
#'
#' Provides smoothing methods for multidimensional scaling-based projections of a given distance matrix.
#'
#' Smoothing is performed using Duchon splines (see \code{\link{Duchon.spline}} for more information).
#'
#' @aliases Predict.matrix.dm.smooth
#' @method Predict.matrix dm.smooth
#' @export
#' @import mgcv
#' @useDynLib msg wood_path
#'
#' @param object a smooth specification object, usually generated by a term \code{s(...,bs="ds",...)}. Note that \code{xt} object is needed, see Details.
#' @param data a list containing just the data (including any \code{by} variable) required by this term, with names corresponding to \code{object$term} (and \code{object$by}). The \code{by} variable is the last element.
Predict.matrix.dm.smooth <- function(object, data){


  mds_obj <- object$xt$mds
  mds.dim <- object$xt$mds.dim

  # use the distance function
  if(is.null(object$xt$dist_fn) | is.null(object$xt$grid)){
    stop("Can't predict without a defined distance function")
  }

  dist_fn <- object$xt$dist_fn
  grid_points <- object$xt$grid
  mds_pred <- insert.mds.generic(mds_obj, data, grid_points, dist_fn)

  mds_pred <- as.data.frame(mds_pred)

  # make some variable names up
  mds.names <- paste("mds-",1:dim(mds_pred)[2],sep="")
  # remove any already in the data
  names(mds_pred) <- mds.names
  object$term <- mds.names
  object$dim <- mds.dim

  Predict.matrix.duchon.spline(object, mds_pred)
}