[adegenet-commits] r375 - pkg/R

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Fri Jun 5 20:02:32 CEST 2009


Author: jombart
Date: 2009-06-05 20:02:31 +0200 (Fri, 05 Jun 2009)
New Revision: 375

Modified:
   pkg/R/haploSim.R
Log:
haplosim now returns ances by label, not by index; labels do not change, while indices become wrong asa we subset the object


Modified: pkg/R/haploSim.R
===================================================================
--- pkg/R/haploSim.R	2009-06-05 14:20:24 UTC (rev 374)
+++ pkg/R/haploSim.R	2009-06-05 18:02:31 UTC (rev 375)
@@ -4,7 +4,7 @@
 ##
 ## N: number of sequences to simulate
 ## mu: mutation rate per nucleotid per year
-## Tmax: periode of time to simulatet
+## Tmax: periode of time to simulate
 ## mean.gen.time, sd.gen.time: average time for transmission and its standard deviation (normal dist)
 ## mean.repro, sd.repro: average number of transmissions and its standard deviation (normal dist)
 ##
@@ -206,12 +206,8 @@
 
 
         ## SHAPE AND RETURN OUTPUT ##
-        nbAncesNAOk <- sum(is.na(res$ances))
-        res$ances <- match(res$ances, rownames(res$seq))
-        if(sum(is.na(res$ances))> nbAncesNAOk){ # in case non-NA ancestors are not in res$seq
-            warning("NA introduced when converting ances to indices, likely indicating a bug")
-        }
-
+        res$ances <- as.character(res$ances)
+        names(res$dates) <- rownames(res$seq)
         class(res) <- "haploSim"
         return(res)
 
@@ -246,11 +242,8 @@
 
 
         ## SHAPE AND RETURN OUTPUT ##
-        nbAncesNAOk <- sum(is.na(res$ances))
-        res$ances <- match(res$ances, rownames(res$seq))
-        if(sum(is.na(res$ances))> nbAncesNAOk){ # in case non-NA ancestors are not in res$seq
-            warning("NA introduced when converting ances to indices, likely indicating a bug")
-        }
+        res$ances <- as.character(res$ances)
+        names(res$dates) <- rownames(res$seq)
 
         class(res) <- "haploSim"
         res$call <- match.call()
@@ -396,3 +389,37 @@
                                      nstep=nstep, step.size=step.size, mu0=mu0,
                                      seq.length=seq.length, rMissDate=rMissDate, ...)
 } # end optimize.seqTrack.haploSim
+
+
+
+
+
+
+
+################
+## plotHaploSim
+################
+plotHaploSim <- function(x, annot=TRUE, dateRange=NULL, col=NULL, bg="grey", add=FALSE, ...){
+
+    ## SOME CHECKS ##
+    if(class(x)!="haploSim") stop("x is not a haploSim object")
+    if(is.null(x$xy)) stop("x does not contain xy coordinates; try to simulate date")
+
+
+    ## CONVERSION TO A SEQTRACK-LIKE OBJECT ##
+    xy <- x$xy
+    res <- list()
+    res$id <- labels(x)
+    res <- as.data.frame(res)
+    res$ances <- x$ances
+    res$weight <- 1 # ??? have to recompute that...
+    res$date <- x$dates
+    res$ances.date <- x$dates[x$ances]
+
+
+    ## CALL TO PLOTSEQTRACK ##
+    out <- plotSeqTrack(res, annot=annot, dateRange=dateRange, col=col, bg=bg, add=add, ...)
+
+    return(invisible(out))
+
+} # end plotHaploSim



More information about the adegenet-commits mailing list