This function fits a point process model
  to an observed point pattern.
  The model may include
  spatial trend, interpoint interaction, and dependence on covariates.
 
  Model-fitting is currently performed by
  the method of maximum pseudolikelihood (Besag, 1975).
  Other options will be added in future versions.
  
  The algorithm is an implementation of the method of
  Baddeley and Turner (2000), which approximates the pseudolikelihood
  by a special type of quadrature sum generalising the Berman-Turner (1992) 
  approximation. 
 
  The argument Q should be either a point pattern
  or a quadrature scheme. If it is a point pattern, it is converted
  into a quadrature scheme.  A quadrature scheme is an object of class "quad"
  (see quad.object)
  which specifies both the data point pattern and the dummy points
  for the quadrature scheme, as well as the quadrature weights
  associated with these points.
  If Q is simply a point pattern
  (of class "ppp", see ppp.object)
  then it is interpreted as specifying the
  data points only; a set of dummy points specified
  by default.dummy() is added, and the default weighting rule is
  invoked to compute the quadrature weights.
 
  The usage of ppm() is closely analogous to the Splus/Rfunctions
  glm and gam.
  The analogy is:
  ll{
  glm 	ppm 
formula 	trend 
family 	interaction
  }
  The point process model to be fitted is specified by the 
  arguments trend and interaction
  which are respectively analogous to
  the formula and family arguments of glm(). 
 
  Systematic effects (spatial trend and/or dependence on 
  spatial covariates) are specified by the argument
  trend. This is an Splus/Rformula object, which may be expressed
  in terms of the Cartesian coordinates x, y,
  the marks marks,
  or the variables in covariates (if supplied), or both.
  It specifies the logarithm of the first order potential
  of the process.
  The formula should not use any names beginning with .mpl
  as these are reserved for internal use.
  If trend is absent or equal to the default, ~1, then
  the model to be fitted is stationary (or at least, its first order 
  potential is constant). 
 
  Stochastic interactions between random points of the point process
  are defined by the argument interaction. This is an object of
  class "interact" which is initialised in a very similar way to the
  usage of family objects in glm and gam.
  See the examples below.
 
  If interaction is missing or NULL, then the model to be fitted
  has no interpoint interactions, that is, it is a Poisson process
  (stationary or nonstationary according to trend). In this case
  the method of maximum pseudolikelihood
  coincides with maximum likelihood.
  To fit a model involving spatial covariates
  other than the Cartesian coordinates $x$ and $y$,
  the values of the covariates should be supplied in the
  argument covariates. 
  Note that it is not sufficient to have observed
  the covariate only at the points of the data point pattern; 
  the covariate must also have been observed at other 
  locations in the window.
  The argument covariates, if supplied, should be either
  a data frame or a list of images.
  If it is a data frame, it must have
  as many rows as there are points in Q.
  The $i$th row of covariates should contain the values of
  spatial variables which have been observed
  at the $i$th point of Q. In this case
  the argument Q must be a quadrature scheme,
  not merely a point pattern.
  Alternatively, covariates may be a list of images,
  each image giving the values of a spatial covariate at 
  a fine grid of locations. 
  The argument covariates should then be a named list,
  with the names corresponding to
  the names of covariates in the model formula trend.
  Each entry in the list must be an image object (of class "im",
  see im.object).
  The software will look up
  the pixel values of each image at the quadrature points.
  In this case the argument Q may be either a quadrature
  scheme or a point pattern.
  The variable names x, y and marks
  are reserved for the Cartesian 
  coordinates and the mark values,
  and these should not be used for variables in covariates.
  
  The argument correction is the name of an edge correction method.
  The default correction="border" specifies the border correction,
  in which the quadrature window (the domain of integration of the 
  pseudolikelihood) is obtained by trimming off a margin of width rbord
  from the observation window of the data pattern.
  Not all edge corrections are implemented (or implementable)
  for arbitrary windows.
  Other options depend on the argument interaction, but these generally 
  include "periodic" (the periodic or toroidal edge correction
  in which opposite edges of a rectangular window are identified)
  and "translate"
  (the translation correction, see Baddeley 1998 and Baddeley and Turner
  2000).
  For pairwise interaction there is also Ripley's isotropic correction,
  identified by "isotropic" or "Ripley".
 
  The fitted point process model returned by this function can be printed 
  (by the print method print.ppm)
  to inspect the fitted parameter values.
  If a nonparametric spatial trend was fitted, this can be extracted using
  the predict method predict.ppm.
  This algorithm approximates the log pseudolikelihood
  by a sum over a finite set of quadrature points.
  Finer quadrature schemes (i.e. those with more
  quadrature points) generally yield a better approximation, at the
  expense of higher computational load. 
  Complete control over the quadrature scheme is possible.
  See quadscheme for an overview.
  Note that the method of maximum pseudolikelihood is
  believed to be inefficient and biased for point processes with strong
  interpoint interactions. In such cases, it is advisable to use
  iterative maximum likelihood
  methods such as Monte Carlo Maximum Likelihood (Geyer, 1999) 
  provided the appropriate simulation algorithm exists. 
  The maximum pseudolikelihood parameter estimate often serves
  as a good initial starting point for these iterative methods. 
  Maximum pseudolikelihood may also be used profitably for
  model selection in the initial phases of modelling.
  The warning message ``Algorithm did not converge in: glm.fit( ...''
  is sometimes obtained. It means that the glm fitting algorithm
  did not converge in a reasonable number of iterations. In most cases
  this can be fixed by taking a finer quadrature scheme.
  See quadscheme.