7.2 R version

The R version uses R's nlm function to maximise the likelihood numerically. nlm does not take box constaints. The idea is therefore to work with parameters transformed onto the whole real line, for example by logit transformation of probabilities (or proportions) (such as $u$ and $m_0$ of the steppingstone model). The user defined function should therefore be able to handle transformed paremeter when the additional argument trans=T. In addition, when the other additional argument convert=T, the user defined function defining the migration pattern should return, instead of a migration matrix, the parameter vector converted back to either the untransformed or transformed state, depending on whether trans=T or trans=F, respectively. Also, when the optional argument start=T, a vector of sensible (untransformed) starting parameter values should be returned. You are free to implement this in any way you like inside your own function. The R version of steppingstone should give you an idea of how to do it:

   
steppingstone <- function(theta,n,trans=F,start=F,convert=F)
{
   if (trans) {theta <- c(ilogit(theta[1],lower=.0001),ilogit(theta[2]))}
   if (convert) {
    if (!trans) {theta <- c(logit(theta[1],lower=.0001),logit(theta[2]))}
    return(theta)
   }   
   if (start) return(c(.1,.1))
   
   M <- matrix(0,ncol=n,nrow=n)
   diag(M) <- 1 - theta[2]
   M[abs(row(M)-col(M))==1] <- theta[2]/2
   M[1,1] <- M[n,n]  <- 1 - theta[2]/2
   M <- (1-theta[1])*M
   return(M)
}

Logit and log transformations and their inverses are provided by the functions logit, ilogit, log, and ilog which are part of migrlib.r. In some cases is necessary to restict the lower bound of the parameter range to some small postive value, e.g. 0.0001 instead of 0 (see Tufto et al., 1998, p. 1979 for details). This lower bound can be passed to logit, ilogit, log, and ilog through the optional lower argument.