[adegenet-commits] r803 - / pkg/R
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Thu Feb 10 16:14:06 CET 2011
Author: jombart
Date: 2011-02-10 16:14:06 +0100 (Thu, 10 Feb 2011)
New Revision: 803
Removed:
adegenet_1.1-0.tar.gz
adegenet_1.1-1.tar.gz
adegenet_1.1-2.tar.gz
Modified:
pkg/R/glFunctions.R
pkg/R/zzz.R
Log:
removed useless tgz; new version glPca.
Deleted: adegenet_1.1-0.tar.gz
===================================================================
(Binary files differ)
Deleted: adegenet_1.1-1.tar.gz
===================================================================
(Binary files differ)
Deleted: adegenet_1.1-2.tar.gz
===================================================================
(Binary files differ)
Modified: pkg/R/glFunctions.R
===================================================================
--- pkg/R/glFunctions.R 2011-02-10 11:35:09 UTC (rev 802)
+++ pkg/R/glFunctions.R 2011-02-10 15:14:06 UTC (rev 803)
@@ -204,104 +204,109 @@
##
## PCA for genlight objects
##
-glPca <- function(x, center=TRUE, scale=FALSE, nf=NULL, loadings=TRUE,
+glPca <- function(x, center=TRUE, scale=FALSE, nf=NULL, loadings=TRUE, alleleAsUnit=FALSE,
useC=TRUE, multicore=require("multicore"), n.cores=NULL){
if(!inherits(x, "genlight")) stop("x is not a genlight object")
- if(multicore && !require(multicore)) stop("multicore package requested but not installed")
- if(multicore && is.null(n.cores)){
- n.cores <- multicore:::detectCores()
- }
+ if(!useC){
+ if(multicore && !require(multicore)) stop("multicore package requested but not installed")
+ if(multicore && is.null(n.cores)){
+ n.cores <- multicore:::detectCores()
+ }
- ## COMPUTE MEANS AND VARIANCES ##
- if(center) {
- vecMeans <- glMean(x, alleleAsUnit=FALSE)
- if(any(is.na(vecMeans))) stop("NAs detected in the vector of means")
- }
- if(scale){
- vecVar <- glVar(x, alleleAsUnit=FALSE)
- if(any(is.na(vecVar))) stop("NAs detected in the vector of variances")
- }
+ ## COMPUTE MEANS AND VARIANCES ##
+ if(center) {
+ vecMeans <- glMean(x, alleleAsUnit=FALSE)
+ if(any(is.na(vecMeans))) stop("NAs detected in the vector of means")
+ }
+ if(scale){
+ vecVar <- glVar(x, alleleAsUnit=FALSE)
+ if(any(is.na(vecVar))) stop("NAs detected in the vector of variances")
+ }
- ## COMPUTE DOT PRODUCTS BETWEEN GENOTYPES ##
- ## to be fast, a particular function is defined for each case of centring/scaling
- myPloidy <- ploidy(x)
+ ## COMPUTE DOT PRODUCTS BETWEEN GENOTYPES ##
+ ## to be fast, a particular function is defined for each case of centring/scaling
- ## NO CENTRING / NO SCALING
- if(!center & !scale){
- dotProd <- function(a,b, ploid.a, ploid.b){ # a and b are two SNPbin objects
- a <- as.integer(a) / ploid.a
- a[is.na(a)] <- 0
- b <- as.integer(b) / ploid.b
- b[is.na(b)] <- 0
- return(sum( a*b, na.rm=TRUE))
+ myPloidy <- ploidy(x)
+
+ ## NO CENTRING / NO SCALING
+ if(!center & !scale){
+ dotProd <- function(a,b, ploid.a, ploid.b){ # a and b are two SNPbin objects
+ a <- as.integer(a) / ploid.a
+ a[is.na(a)] <- 0
+ b <- as.integer(b) / ploid.b
+ b[is.na(b)] <- 0
+ return(sum( a*b, na.rm=TRUE))
+ }
}
- }
- ## CENTRING / NO SCALING
- if(center & !scale){
- dotProd <- function(a,b, ploid.a, ploid.b){ # a and b are two SNPbin objects
- a <- as.integer(a) / ploid.a
- a[is.na(a)] <- vecMeans[is.na(a)]
- b <- as.integer(b) / ploid.b
- b[is.na(b)] <- vecMeans[is.na(b)]
- return(sum( (a-vecMeans) * (b-vecMeans), na.rm=TRUE) )
+ ## CENTRING / NO SCALING
+ if(center & !scale){
+ dotProd <- function(a,b, ploid.a, ploid.b){ # a and b are two SNPbin objects
+ a <- as.integer(a) / ploid.a
+ a[is.na(a)] <- vecMeans[is.na(a)]
+ b <- as.integer(b) / ploid.b
+ b[is.na(b)] <- vecMeans[is.na(b)]
+ return(sum( (a-vecMeans) * (b-vecMeans), na.rm=TRUE) )
+ }
}
- }
- ## NO CENTRING / SCALING (odd option...)
- if(!center & scale){
- dotProd <- function(a,b, ploid.a, ploid.b){ # a and b are two SNPbin objects
- a <- as.integer(a) / ploid.a
- a[is.na(a)] <- 0
- b <- as.integer(b) / ploid.b
- b[is.na(b)] <- 0
- return(sum( (a*b)/vecVar, na.rm=TRUE))
+ ## NO CENTRING / SCALING (odd option...)
+ if(!center & scale){
+ dotProd <- function(a,b, ploid.a, ploid.b){ # a and b are two SNPbin objects
+ a <- as.integer(a) / ploid.a
+ a[is.na(a)] <- 0
+ b <- as.integer(b) / ploid.b
+ b[is.na(b)] <- 0
+ return(sum( (a*b)/vecVar, na.rm=TRUE))
+ }
}
- }
- ## CENTRING / SCALING
- if(center & scale){
- dotProd <- function(a,b, ploid.a, ploid.b){ # a and b are two SNPbin objects
- a <- as.integer(a) / ploid.a
- a[is.na(a)] <- vecMeans[is.na(a)]
- b <- as.integer(b) / ploid.b
- a[is.na(a)] <- vecMeans[is.na(a)]
- return( sum( ((a-vecMeans)*(b-vecMeans))/vecVar, na.rm=TRUE ) )
+ ## CENTRING / SCALING
+ if(center & scale){
+ dotProd <- function(a,b, ploid.a, ploid.b){ # a and b are two SNPbin objects
+ a <- as.integer(a) / ploid.a
+ a[is.na(a)] <- vecMeans[is.na(a)]
+ b <- as.integer(b) / ploid.b
+ a[is.na(a)] <- vecMeans[is.na(a)]
+ return( sum( ((a-vecMeans)*(b-vecMeans))/vecVar, na.rm=TRUE ) )
+ }
}
- }
- ## COMPUTE ALL POSSIBLE DOT PRODUCTS (XX^T / n) ##
- allComb <- combn(1:nInd(x), 2)
- if(multicore){
- allProd <- unlist(mclapply(1:ncol(allComb), function(i) dotProd(x at gen[[allComb[1,i]]], x at gen[[allComb[2,i]]], myPloidy[allComb[1,i]], myPloidy[allComb[2,i]]),
- mc.cores=n.cores, mc.silent=TRUE, mc.cleanup=TRUE, mc.preschedule=FALSE))
- } else {
- allProd <- unlist(lapply(1:ncol(allComb), function(i) dotProd(x at gen[[allComb[1,i]]], x at gen[[allComb[2,i]]], myPloidy[allComb[1,i]], myPloidy[allComb[2,i]]) ))
- }
- allProd <- allProd / nInd(x) # assume uniform weights
+ ## COMPUTE ALL POSSIBLE DOT PRODUCTS (XX^T / n) ##
+ allComb <- combn(1:nInd(x), 2)
+ if(multicore){
+ allProd <- unlist(mclapply(1:ncol(allComb), function(i) dotProd(x at gen[[allComb[1,i]]], x at gen[[allComb[2,i]]], myPloidy[allComb[1,i]], myPloidy[allComb[2,i]]),
+ mc.cores=n.cores, mc.silent=TRUE, mc.cleanup=TRUE, mc.preschedule=FALSE))
+ } else {
+ allProd <- unlist(lapply(1:ncol(allComb), function(i) dotProd(x at gen[[allComb[1,i]]], x at gen[[allComb[2,i]]], myPloidy[allComb[1,i]], myPloidy[allComb[2,i]]) ))
+ }
+ allProd <- allProd / nInd(x) # assume uniform weights
- ## shape result as a matrix
- attr(allProd,"Size") <- nInd(x)
- attr(allProd,"Diag") <- FALSE
- attr(allProd,"Upper") <- FALSE
- class(allProd) <- "dist"
- allProd <- as.matrix(allProd)
+ ## shape result as a matrix
+ attr(allProd,"Size") <- nInd(x)
+ attr(allProd,"Diag") <- FALSE
+ attr(allProd,"Upper") <- FALSE
+ class(allProd) <- "dist"
+ allProd <- as.matrix(allProd)
- ## compute the diagonal
- if(multicore){
- temp <- unlist(mclapply(1:nInd(x), function(i) dotProd(x at gen[[i]], x at gen[[i]], myPloidy[i], myPloidy[i]),
- mc.cores=n.cores, mc.silent=TRUE, mc.cleanup=TRUE, mc.preschedule=FALSE))/nInd(x)
- } else {
- temp <- unlist(lapply(1:nInd(x), function(i) dotProd(x at gen[[i]], x at gen[[i]], myPloidy[i], myPloidy[i]) ))/nInd(x)
+ ## compute the diagonal
+ if(multicore){
+ temp <- unlist(mclapply(1:nInd(x), function(i) dotProd(x at gen[[i]], x at gen[[i]], myPloidy[i], myPloidy[i]),
+ mc.cores=n.cores, mc.silent=TRUE, mc.cleanup=TRUE, mc.preschedule=FALSE))/nInd(x)
+ } else {
+ temp <- unlist(lapply(1:nInd(x), function(i) dotProd(x at gen[[i]], x at gen[[i]], myPloidy[i], myPloidy[i]) ))/nInd(x)
+ }
+ diag(allProd) <- temp
+ } else { # === use C computations ====
+ glDotProd(x, center=center, scale=scale, alleleAsUnit=alleleAsUnit)
}
- diag(allProd) <- temp
## PERFORM THE ANALYSIS ##
Modified: pkg/R/zzz.R
===================================================================
--- pkg/R/zzz.R 2011-02-10 11:35:09 UTC (rev 802)
+++ pkg/R/zzz.R 2011-02-10 15:14:06 UTC (rev 803)
@@ -2,8 +2,7 @@
#.initAdegenetClasses()
#.initAdegenetUtils()
library.dynam("adegenet", pkg, lib)
- startup.txt <- " ==========================\n adegenet 1.2-9 is loaded\n ==========================\n\n - to start, type '?adegenet'\n - to browse adegenet website, type 'adegenetWeb()'\n - to post questions/comments: adegenet-forum at lists.r-forge.r-project.org\n\n"
+ startup.txt <- " ==========================\n adegenet 1.3-0 is loaded\n ==========================\n\n - to start, type '?adegenet'\n - to browse adegenet website, type 'adegenetWeb()'\n - to post questions/comments: adegenet-forum at lists.r-forge.r-project.org\n\n"
packageStartupMessage(startup.txt)
-
}
More information about the adegenet-commits
mailing list