[Mattice-commits] r59 - in pkg: R man

[noreply at r-forge.r-project.org]

Author: andrew_hipp
Date: 2008-12-01 17:37:00 +0100 (Mon, 01 Dec 2008)
New Revision: 59

Modified:
   pkg/R/batchHansen.R
   pkg/man/runBatchHansen.Rd
Log:
Updates to allow user to save 'hansen' objects called by runBatchHansen

Modified: pkg/R/batchHansen.R
===================================================================
--- pkg/R/batchHansen.R	2008-11-25 19:10:48 UTC (rev 58)
+++ pkg/R/batchHansen.R	2008-12-01 16:37:00 UTC (rev 59)
@@ -2,8 +2,8 @@
 # FUNCTIONS FOR PERFORMING A SERIES OF OU ANALYSES ON A BATCH OF TREES
 # ---------------------------------------------------------------------
 ## Changes needed:
-## 2. measurement error portions need to be fixed
-## 5. In a better world, allow graphical selection of subtrees to test on a single tree, then extract defining taxa
+## - measurement error portions need to be fixed
+## - In a better world, allow graphical selection of subtrees to test on a single tree, then extract defining taxa
 ##    based on those nodes, using locator() or something like it.
 
 runBatchHansen <-
@@ -17,7 +17,7 @@
 #  "cladeMembersList" = list of vectors containing names of the members of each clade (except for the root of the tree)
 #  "brown" = whether to analyse the data under a Brownian motion model
 #  "..." = additional arguments to pass along to hansen
-function(ouchTrees, characterStates, cladeMembersList, nodeNames = NULL, maxNodes = NULL, regimeTitles = NULL, brown = F, rescale = 1, ...) {
+function(ouchTrees, characterStates, cladeMembersList, filePrefix = NULL, di = NULL, nodeNames = NULL, maxNodes = NULL, regimeTitles = NULL, brown = F, rescale = 1, ...) {
   ## do all the objects in ouchTrees inherit ouchtree?
   if(is(ouchTrees,'ouchtree')) ouchTrees <- list(ouchTrees)
   treeCheck <- unlist(lapply(ouchTrees, function(x) is(x,'ouchtree')))
@@ -47,10 +47,12 @@
     else dataFlag <- 'named'
     if(stopFlag) stop("Correct discrepancies between trees and data and try again!")
     }
-
+  if(!identical(di, NULL)) dir.create(di)
   ar = regimeVectors(ouchTrees, cladeMembersList, maxNodes)
   hansenBatch <- list(length(ouchTrees))
   for (i in 1:length(ouchTrees)) {
+    if(!identical(filePrefix, NULL)) fP <- paste(filePrefix, ".t", i, ".", sep = "")
+    if(!identical(di, NULL)) fP <- paste(di, "/", fP, sep = "")
     tree <- ouchTrees[[i]]
     if(identical(regimeTitles, NULL)) {
       regimeTitles <- as.character(1:length(ar$regList[[i]]))
@@ -71,7 +73,7 @@
     ## send it off to batchHansen and just stick the results in hansenBatch... this won't work as the number of regimes gets large, 
     ##   so there should be some option here to just hang onto the coefficients for each run (i.e., hang onto 'coef(hansen(...))' rather than 'hansen(...)')
     ##   there could also be an option to save the entire object as a series of files in addition to hanging onto 
-    hansenBatch[[i]] <- batchHansen(tree, dataIn, ar$regList[[i]], regimeTitles, brown, ...)
+    hansenBatch[[i]] <- batchHansen(tree, dataIn, ar$regList[[i]], regimeTitles, brown, fP, ...)
     message(paste("Tree",i,"of",length(ouchTrees),"complete", "\n-----------------------------"))
   }
     
@@ -87,7 +89,7 @@
 #  "regimesList" = list of regime-paintings as output from regimeVectors
 #  "scalingFactor" = factor to multiply against (times / max(times)) -- choose based on trial analyses
 # Value: a matrix with nrow = regimes (+ 1 if brownian model is included) and columns for u, d.f., all estimated parameters, LRvsBM, AIC, and AIC weight
-function(tree, data, regimesList, regimeTitles, brown, ...) {
+function(tree, data, regimesList, regimeTitles, brown, filePrefix, ...) {
   n <- tree at nterm
   ## set up a matrix that returns lnL, K, sigmasq, theta0, and alpha for every model; thetas will go along into a list that is indexed by model
   hansenOptima <- list(length(regimeTitles))
@@ -109,6 +111,7 @@
       message(paste("Running regime",i))
       ## at this point, the user has to give an initial alpha and sigma for hansen to search on... this should be relaxed
       ha = hansen(data, tree, regimesList[[i]], ...)
+      if(!identical(filePrefix, NULL)) save(ha, paste(filePrefix, 'r', i, '.Rdata', sep = ""))
       treeData[i, ] <- unlist(summary(ha)[haVars])
       hansenOptima[[i]] <- summary(ha)$optima[[1]]
       }

Modified: pkg/man/runBatchHansen.Rd
===================================================================
--- pkg/man/runBatchHansen.Rd	2008-11-25 19:10:48 UTC (rev 58)
+++ pkg/man/runBatchHansen.Rd	2008-12-01 16:37:00 UTC (rev 59)
@@ -1,13 +1,12 @@
 \name{runBatchHansen}
 \alias{runBatchHansen}
-\title{Analyze Ornstein-Uhlenbeck models over a set of trees, with changes in regime corresponding to permutations
-       of the nodes fed to it.}
+\title{Analyze Ornstein-Uhlenbeck models over a set of trees}
 \description{
   \code{runBatchHansen} simulates the evolution of a single character on an \code{ape}-style user tree under arbitrarily
   complex Ornstein-Uhlenbeck models.
 }
 \usage{
-runBatchHansen(ouchTrees, characterStates, cladeMembersList, nodeNames = NULL, maxNodes = NULL, regimeTitles = NULL, brown = F, rescale = 1, ...)}
+runBatchHansen(ouchTrees, characterStates, cladeMembersList, filePrefix = NULL, di = NULL, nodeNames = NULL, maxNodes = NULL, regimeTitles = NULL, brown = F, rescale = 1, ...)}
 \arguments{
   \item{ouchTrees}{
 	A list of \code{ouchtree} style trees, or a single tree
@@ -21,6 +20,13 @@
 	A list of tip-name vectors, each defining a node of interest in the tree. Because these vectors are checked for monophyly on each tree, each
 	node must by defined by all tips descendent from that node.
     }
+  \item{filePrefix}{
+  	A file prefix for saving hansen objects. If left at default (\code{NULL}), no hansen objects will be saved. Typically 
+  	this is fine, but if you want to extract parameters from the hansen objects, or do any bootstrapping, you'll need these.
+    }
+  \item{di}{
+	A directory in which to save hansen objects. Not required, but strongly recommended if you are saving hansen objects.  
+    }
   \item{nodeNames}{
 	A vector of names for the nodes being analyzed. Not implemented. Optional.
     }
@@ -36,19 +42,20 @@
   \item{rescale}
 	Enter a value > 0 to rescale all trees to length = \code{rescale}.
     }
-  \item{...}
+  \item{...}{
     Arguments passed along to \code{hansen}. As written now, you have to include a beginning value for \code{sigma}
     and \code{alpha} for the optimizer to begin with. Other arguments are described in documentation for \code{hansen}.
     To begin, it is easiest to set \code{sigma} = 1, \code{alpha} = 1.
     }
 }
 \details{
-    This function is the primary function for estimating the 
+    This function is the primary function for estimating the probability of a change in with changes in regime at the nodes
+    specified, with changes corresponding to permutations of the nodes fed to the function. See vignette for additional details.
 }
 \value{
   XXX
   }
-\author{Andrew Hipp ahipp at mortonarb.org}
+\author{Andrew Hipp <ahipp at mortonarb.org>}
 \seealso{
    \code{batchHansen} for running analysis on a single tree under an explicit regime,
    \code{summary.hansenBatch} for summarizing data.