[Mattice-commits] r66 - in pkg: . R
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Tue Dec 2 19:27:35 CET 2008
Author: andrew_hipp
Date: 2008-12-02 19:27:35 +0100 (Tue, 02 Dec 2008)
New Revision: 66
Added:
pkg/TODO
Modified:
pkg/R/regimes.R
pkg/R/summarizingAnalyses.R
pkg/R/treeTraversal.R
Log:
misc. cleanup
Modified: pkg/R/regimes.R
===================================================================
--- pkg/R/regimes.R 2008-12-02 14:29:11 UTC (rev 65)
+++ pkg/R/regimes.R 2008-12-02 18:27:35 UTC (rev 66)
@@ -19,12 +19,12 @@
# "nodeMatrix" = matrix of trees (rows) by nodes (columns) indicating what nodes are present on which trees
# 19 nov 08: changing to accept a list of trees and trimmed down greatly
function(ouchTrees, cladeMembersList, maxNodes = NULL) {
-
nnode <- length(cladeMembersList)
regMatrix <- regimeMatrix(n = nnode, maxNodes = maxNodes)
apr = regimeMaker(ouchTrees, regMatrix, cladeMembersList)
outdata <- list(regList = apr$regList, regMatrix = apr$regMatrix, nodeMatrix = apr$nodeMatrix)
- return(outdata) }
+ return(outdata)
+}
paintBranches <-
# Paints branches with regimes changing at nodes specified
Modified: pkg/R/summarizingAnalyses.R
===================================================================
--- pkg/R/summarizingAnalyses.R 2008-12-02 14:29:11 UTC (rev 65)
+++ pkg/R/summarizingAnalyses.R 2008-12-02 18:27:35 UTC (rev 66)
@@ -8,12 +8,7 @@
## items in output: hansens, regimeList, regimeMatrix
## the summary will eventually sum weights over all nodes over all trees
- # Unimplemented summary ideas
- # - Check whether there is a single tree
- # - if so, return everything below, + a model-averaged set of thetas indexed according to branches
- # all that's needed to do now is model average over the theta matrix created for each tree
-
- # 0. Get information criterion weights for all models
+ # Get information criterion weights for all models
icObject <- informationCriterion.hansenBatch(hansenBatch)
nmodels <- dim(hansenBatch$hansens[[1]])[1]
ntrees <- length(hansenBatch$hansens)
@@ -46,6 +41,7 @@
weighted.mean,
w = w,
na.rm = T)
+
## the lines below made the weights on branches ignore the fact that the Brownian motion model was part of the
## model set; however, I've removed them b/c support for the Brownian motion model does (and should) contribute
## to reduced probability of change at any of the nodes. You can uncomment them if you feel differently.
@@ -65,7 +61,7 @@
# acting being a proper probability distribution, albeit one that confounds clade support
# with model support.
- # 1. sum over nodes
+ # sum over nodes
nodes <- dimnames(hansenBatch$regMatrix$overall)[[2]] # grab the overall regMatrix, which includes all possible nodes, no matter which trees do or don't have them
nodeWeightsMatrix.unnormalized <- matrix(NA, nrow = length(matrixRows), ncol = length(nodes), dimnames = list(matrixRows, nodes)) # value: nodeWeightsMatrix
nodeWeightsMatrix.allNodes <- nodeWeightsMatrix.unnormalized # same dimensions, another empty matrix to fill up
@@ -81,7 +77,7 @@
else nodeWeightsMatrix.unnormalized[, nodes[i]] <- apply(modelsMatrixSubset, 2, sum)
}
- # 2. sum over number of parameters
+ # sum over number of parameters
# create a vector of sums that tells us how many categories there are for each model: dof = sum(nodes) + 1 [because a node indicates a change in
# regime, thus the total number of thetas = nodes + 1] + alpha + sigma = sum(nodes) + 3; for Brownian motion model, dof = 2
nodeSums <- apply(hansenBatch$regMatrix$overall, 1, sum) + 3
@@ -108,6 +104,7 @@
}
print.hansenSummary <- function(hansenSummary) {
+## This just formats a hansenSummary object so that it is readable on the screen; you can still store the summary object and extract elements as needed
message(paste("Summarizing hansenBatch analyses over", length(hansenSummary$modelsMatrix), "trees and", dim(hansenSummary$modelsMatrix[[1]])[1], "models"))
message("ESTIMATED SUPPORT FOR CHANGES OCCURRING AT DESIGNATED NODES\nThese estimates are only valid if (1) the maximum number of regimes permitted approximates the actual maximum;\n(2) nodes at which changes actually occurred were included among the nodes being tested; and\n(3) any matrix you may have utilized to conduct analysis was balanced, such that all nodes are present in the same number of models.\n\n")
message("Averaged only over models containing that node:")
Modified: pkg/R/treeTraversal.R
===================================================================
--- pkg/R/treeTraversal.R 2008-12-02 14:29:11 UTC (rev 65)
+++ pkg/R/treeTraversal.R 2008-12-02 18:27:35 UTC (rev 66)
@@ -2,23 +2,6 @@
# FUNCTIONS FOR TRAVERSING AN S4 OUCH TREE #
# ------------------------------------------
-# Modified from functions used in Hipp 2007 Evolution paper
-# Initially written for ouch v 1.2-4
-# 10 November 2008: changed everything to operate on an ouchtree object (ouch >= v2), otherwise functions the same;
-# checked on ouch v2.4-2
-# functions included in this file:
-# 1. paintBranches
-# 2. mrcaOUCH
-# 3. ancestorLine
-# 4. allPossibleRegimes
-# 5. regimeVectors
-
-# 19 nov 08: moved all regime-involved routines to regimes.R
-
-# To do:
-# 2. make allPossibleRegimes more efficient when maxNodes < length(nodes)
-
-
isMonophyletic <- function(tree, taxa) {
# returns T or F on whether a group of taxa is monophyletic in an ouch tree
if(length(taxa) == 1) return(taxa %in% tree at nodelabels[tree at term])
@@ -27,7 +10,6 @@
nodeDescendents <- function(tree, startNode) {
## Recursive function to find all the descendents of a node on an 'ouchtree' object
-## a bit clunky as written
startNode <- as.character(startNode) # just to be safe
daughterBranches <- as.character(tree at nodes[tree at ancestors %in% startNode])
nodeNames <- tree at nodelabels[tree at nodes %in% daughterBranches]
@@ -45,8 +27,7 @@
# Arguments:
# "node" "ancestor" "times" "species" = the standard tree specification vectors of the OUCH-style tree
# "cladeVector" = vector of species for which you want to find the most recent common ancestor
-# Value: the node number (as an integer) of the most recent common ancestor
-# Works! 3-31-06
+# Value: the node number of the most recent common ancestor
function(cladeVector, tree) {
## ------------------ begin ouchtree block -----------------
## check to see if tree inherits 'ouchtree'
@@ -79,31 +60,16 @@
ancestorLine <-
# Creates a vector of ancestral nodes for a tip
# Arguments:
-# CHANGED to "tree" 10 nov 08: "node" and "ancestor" = the standard tree specification vectors of the OUCH-style tree
+# "tree" = an ouch-style (S4) tree
# "tip" = the tip node to trace back
# Value: a vector of nodes leading from a tip to the root
# 10 nov 08: changed to just grab tree at lineages and make a vector that fits the old code
-# this should really be exluded at some point
function(tip, tree) {
- ## ------------------ begin ouchtree block -----------------
## check to see if tree inherits 'ouchtree'
if (!is(tree,'ouchtree'))
stop(paste('This function has been rewritten to use the new S4 ', sQuote('ouchtree'), ' class.',
'\nYou can generate a tree of this class by calling ', sQuote('ouchtree()'), '.', sep = ""))
- ## get the vectors we need:
- #ancestor <- tree at ancestors # class = "character"
- #node <- tree at nodes # class = "character"
- #species <- tree at nodelabels # class = "character" -- note that nodelabels is more general than this indicates and the name should be changed throughout at some point
- #times <- tree at times # class = "numeric"
- ## ------------------ end ouchtree block -------------------
tip <- as.numeric(tip)
- #nodesVector = vector("character")
- #counter = 0
- #repeat {
- # if (is.na(tip)) break
- # counter = counter + 1
- # nodesVector[counter] = ancestor[tip]
- # tip = ancestor[tip] }
nodesVector <- c(as.character(tree at lineages[[tip]][2:length(tree at lineages[[tip]])]), NA)
return(nodesVector)
}
\ No newline at end of file
Added: pkg/TODO
===================================================================
--- pkg/TODO (rev 0)
+++ pkg/TODO 2008-12-02 18:27:35 UTC (rev 66)
@@ -0,0 +1 @@
+make allPossibleRegimes more efficient when maxNodes < length(nodes)
\ No newline at end of file
More information about the Mattice-commits
mailing list