https://github.com/cran/spatstat
Raw File
Tip revision: 3aca716ce2576a0dab83f08052acd47afed8ee6a authored by Adrian Baddeley on 29 February 2012, 00:00:00 UTC
version 1.25-4
Tip revision: 3aca716
interact.R
#
#	interact.S
#
#
#	$Revision: 1.20 $	$Date: 2012/01/18 09:31:39 $
#
#	Class 'interact' representing the interpoint interaction
#               of a point process model
#              (e.g. Strauss process with a given threshold r)
#
#       Class 'isf' representing a generic interaction structure
#              (e.g. pairwise interactions)
#
#	These do NOT specify the "trend" part of the model,
#	only the "interaction" component.
#
#               The analogy is:
#
#                       glm()             ppm()
#
#                       model formula     trend formula
#
#                       family            interaction
#
#               That is, the 'systematic' trend part of a point process
#               model is specified by a 'trend' formula argument to ppm(),
#               and the interpoint interaction is specified as an 'interact'
#               object.
#
#       You only need to know about these classes if you want to
#       implement a new point process model.
#
#       THE DISTINCTION:
#       An object of class 'isf' describes an interaction structure
#       e.g. pairwise interaction, triple interaction,
#       pairwise-with-saturation, Dirichlet interaction.
#       Think of it as determining the "order" of interaction
#       but not the specific interaction potential function.
#
#       An object of class 'interact' completely defines the interpoint
#       interactions in a specific point process model, except for the
#       regular parameters of the interaction, which are to be estimated
#       by ppm() or otherwise. An 'interact' object specifies the values
#       of all the 'nuisance' or 'irregular' parameters. An example
#       is the Strauss process with a given, fixed threshold r
#       but with the parameters beta and gamma undetermined.
#
#       DETAILS:
#
#       An object of class 'isf' contains the following:
#
#	     $name               Name of the interaction structure         
#                                        e.g. "pairwise"
#
#	     $print		 How to 'print()' this object
#				 [A function; invoked by the 'print' method
#                                 'print.isf()']
#
#            $eval               A function which evaluates the canonical
#                                sufficient statistic for an interaction
#                                of this general class (e.g. any pairwise
#                                interaction.)
#
#       If lambda(u,X) denotes the conditional intensity at a point u
#       for the point pattern X, then we assume
#                  log lambda(u, X) = theta . S(u,X)
#       where theta is the vector of regular parameters,
#       and we call S(u,X) the sufficient statistic.
#
#       A typical calling sequence for the $eval function is
#
#            (f$eval)(X, U, E, potentials, potargs, correction)
#
#       where X is the data point pattern, U is the list of points u
#       at which the sufficient statistic S(u,X) is to be evaluated,
#       E is a logical matrix equivalent to (X[i] == U[j]),
#       $potentials defines the specific potential function(s) and
#       $potargs contains any nuisance/irregular parameters of these
#       potentials [the $potargs are passed to the $potentials without
#       needing to be understood by $eval.]
#       $correction is the name of the edge correction method.
#
#
#       An object of class 'interact' contains the following:
#
#
#            $name               Name of the specific potential
#                                        e.g. "Strauss"
#
#            $family              Object of class "isf" describing
#                                the interaction structure
#
#            $pot	         The interaction potential function(s)
#                                -- usually a function or list of functions.
#                                (passed as an argument to $family$eval)
#
#            $par                list of any nuisance/irregular parameters
#                                (passed as an argument to $family$eval)
#
#            $parnames           vector of long names/descriptions
#                                of the parameters in 'par'
#
#            $init()             initialisation action
#                                or NULL indicating none required
#
#            $update()           A function to modify $par
#                                [Invoked by 'update.interact()']
#                                or NULL indicating a default action
#
#	     $print		 How to 'print()' this object
#				 [Invoked by 'print' method 'print.interact()']
#                                or NULL indicating a default action
#
# --------------------------------------------------------------------------

print.isf <- function(x, ...) {
  verifyclass(x, "isf")
  if(!is.null(x$print))
    (x$print)(x)
  invisible(NULL)
}

print.interact <- function(x, ..., family=TRUE, brief=FALSE) {
  verifyclass(x, "interact")
  if(!is.null(x$print))
     (x$print)(x)
  else {
    # default
    if(family && !brief) 
      print.isf(x$family)
    if(!brief)
      cat("Interaction:")
    cat(paste(x$name, "\n"))
    # just print the parameter names and their values
    cat(paste(x$parnames, ":\t", x$par, "\n", sep=""))
  }
  invisible(NULL)
}

is.interact <- function(x) { inherits(x, "interact") }

update.interact <- function(object, ...) {
  verifyclass(object, "interact")
  if(!is.null(object$update))
    (object$update)(object, ...)
  else {
    # Default
    # First update the version
    if(outdated.interact(object))
      object <- reincarnate.interact(object)
    # just match the arguments in "..."
    # with those in object$par and update them
    want <- list(...)
    if(length(want) > 0) {
      m <- match(names(want),names(object$par))
      nbg <- is.na(m)
      if(any(nbg)) {
        which <- paste((names(want))[nbg])
        warning(paste("Arguments not matched: ", which))
      }
      m <- m[!nbg]
      object$par[m] <- want
    }
    # call object's own initialisation routine
    if(!is.null(object$init))
      (object$init)(object)
    object
  }    
}

  
is.poisson.interact <- function(x) {
  verifyclass(x, "interact")
  is.null(x$family)
}


# Test whether interact object was made by an older version of spatstat

outdated.interact <- function(object) {
  ver <- object$version
  older <- is.null(ver) || (package_version(ver) < versionstring.spatstat())
  return(older)
}


# Test whether the functions in the interaction object
# expect the coefficient vector to contain ALL coefficients,
# or only the interaction coefficients.
# This change was introduced in 1.11-0, at the same time
# as interact objects were given version numbers.

newstyle.coeff.handling <- function(object) {
  stopifnot(inherits(object, "interact"))  
  ver <- object$version
  old <- is.null(ver) || (package_version(ver) < "1.11")
  return(!old)
}

# ######
#
# Re-create an interact object in the current version of spatstat
#
# 

reincarnate.interact <- function(object) {
  # re-creates an interact object in the current version of spatstat

  if(!is.null(object$update)) {
    newobject <- (object$update)(object)
    return(newobject)
  }
  
  par <- object$par
  pot <- object$pot
  name <- object$name
  
  # get creator function
  creator <- object$creator
  if(is.null(creator)) {
    # old version: look up list
    creator <- .Spatstat.Old.InteractionList[[name]]
    if(is.null(creator))
      stop(paste("Don't know how to update", sQuote(name),
                 "to current version of spatstat"))
  }
  if(is.character(creator))
    creator <- get(creator)
  if(!is.function(creator) && !is.expression(creator))
    stop("Internal error: creator is not a function or expression")

  # call creator
  if(is.expression(creator)) 
    newobject <- eval(creator)
  else {
    # creator is a function
  
    # It's assumed that the creator function's arguments are
    # either identical to components of 'par' (the usual case)
    # or to one of the components of the object itself (Ord, Saturated)
    # or to printfun=object$print (Pairwise).
  
    argnames <- names(formals(creator))
    available <- append(par, object)
    available <- append(available, list(printfun=object$print))
    ok <- argnames %in% names(available)
    if(!all(ok))
      stop(paste("Internal error:",
                 ngettext(sum(!ok), "argument", "arguments"),
                 paste(sQuote(argnames[!ok]), collapse=", "),
                 "in creator function were not understood"))
    newobject <- do.call(creator, available[argnames])
  }
  
  if(!inherits(newobject, "interact"))
    stop("Internal error: creator did not return an object of class interact")

  return(newobject)
}


# This list is necessary to deal with older formats of 'interact' objects
# which did not include the creator name

.Spatstat.Old.InteractionList <-
  list("Diggle-Gratton process"    = "DiggleGratton",
       "Geyer saturation process"  = "Geyer",
       "Lennard-Jones potential"   = "LennardJones",
       "Multitype Strauss process" = "MultiStrauss",
       "Multitype Strauss Hardcore process" = "MultiStraussHard",
       "Ord process with threshold potential"="OrdThresh",
       "Piecewise constant pairwise interaction process"="PairPiece",
       "Poisson process"           = "Poisson",
       "Strauss process"           = "Strauss",
       "Strauss - hard core process" = "StraussHard",
       "Soft core process" = "Softcore",
       # weird ones:
       "Ord process with user-defined potential" = expression(Ord(object$pot)),
       "Saturated process with user-defined potential"
          =expression(Saturated(object$pot)),
       "user-defined pairwise interaction process"=
       expression(
           Pairwise(object$pot,
                    par=object$par,
                    parnames=object$parnames,
                    printfun=object$print))
           
     )
       
as.interact <- function(object) {
  UseMethod("as.interact")
}

as.interact.interact <- function(object) {
  verifyclass(object, "interact")
  return(object)
}


#### internal code for streamlining initialisation of interactions
#
#    x should be a partially-completed 'interact' object
#

instantiate.interact <- function(x, par) {
  if(is.character(x$family)) x$family <- get(x$family)
  # set parameter values
  x$par    <- par
  # validate parameters
  x$init(x)
  x$version <- versionstring.spatstat()
  return(x)
}
back to top