[adegenet-commits] r849 - in pkg: . R
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Thu Mar 10 17:42:42 CET 2011
Author: jombart
Date: 2011-03-10 17:42:42 +0100 (Thu, 10 Mar 2011)
New Revision: 849
Added:
pkg/R/glSim.R
Modified:
pkg/DESCRIPTION
pkg/R/glHandle.R
Log:
put glSim in its own source file.
Modified: pkg/DESCRIPTION
===================================================================
--- pkg/DESCRIPTION 2011-03-10 16:37:30 UTC (rev 848)
+++ pkg/DESCRIPTION 2011-03-10 16:42:42 UTC (rev 849)
@@ -10,6 +10,6 @@
Suggests: genetics, hierfstat, spdep, tripack, ape, pegas, graph, RBGL, seqinr, multicore
Depends: methods, MASS, ade4
Description: Classes and functions for genetic data analysis within the multivariate framework.
-Collate: classes.R basicMethods.R handling.R auxil.R setAs.R SNPbin.R glHandle.R glFunctions.R find.clust.R hybridize.R scale.R fstat.R import.R seqTrack.R chooseCN.R genind2genpop.R loadingplot.R sequences.R gstat.randtest.R makefreq.R colorplot.R monmonier.R spca.R coords.monmonier.R haploGen.R old2new.R spca.rtests.R dapc.R haploPop.R PCtest.R dist.genpop.R Hs.R propShared.R export.R HWE.R propTyped.R inbreeding.R glPlot.R zzz.R
+Collate: classes.R basicMethods.R handling.R auxil.R setAs.R SNPbin.R glHandle.R glFunctions.R glSim.R find.clust.R hybridize.R scale.R fstat.R import.R seqTrack.R chooseCN.R genind2genpop.R loadingplot.R sequences.R gstat.randtest.R makefreq.R colorplot.R monmonier.R spca.R coords.monmonier.R haploGen.R old2new.R spca.rtests.R dapc.R haploPop.R PCtest.R dist.genpop.R Hs.R propShared.R export.R HWE.R propTyped.R inbreeding.R glPlot.R zzz.R
License: GPL (>=2)
LazyLoad: yes
Modified: pkg/R/glHandle.R
===================================================================
--- pkg/R/glHandle.R 2011-03-10 16:37:30 UTC (rev 848)
+++ pkg/R/glHandle.R 2011-03-10 16:42:42 UTC (rev 849)
@@ -235,106 +235,8 @@
-##########
-## glSim
-##########
-glSim <- function(n.ind, n.snp.nonstruc, n.snp.struc=0, grp.size=round(n.ind/2), ploidy=1, alpha=0,
- block.size=NULL){
- ## BASIC CHECKS ##
- if( any(c(n.ind, n.snp.nonstruc+n.snp.struc) <1)) stop("null numbers of individuals and/or SNPs requested")
- ## alpha parameter
- if(alpha>0.5){
- alpha <- 0.5
- warning("alpha cannot exceed 0.5 - changing the value to 0.5 (total forced asymmetry)")
- }
- if(alpha<0){
- alpha <- 0
- warning("alpha cannot be lower than 0 - changing the value to 0 (no forced asymmetry)")
- }
- ## handle block size
- if(is.null(block.size)){
- block.size <- n.snp.nonstruc + n.snp.struc
- }
-
- ## handle group sizes
- grpA.size <- grp.size
- if(grpA.size >= n.ind) stop("grpA.size is >= n.ind")
- grpB.size <- n.ind - grpA.size
-
-
- ## AUXIL FUNCTIONS ##
- ## draw p snp for i indiv and convert into a genlight - no structure
- f1 <- function(n,p){
- temp <- sapply(1:p, function(i) rbinom(n, ploidy, runif(1)))
- if(n==1) {temp <- matrix(temp,nrow=1)}
- return(new("genlight", temp, ploidy=ploidy))
- }
-
- ## draw p snp for i indiv and convert into a genlight - differences between 2 groups
- if(n.snp.struc > 0){
- f2 <- function(n,p){
- probA <- runif(p, min=0, max=0.5-alpha)
- probB <- 1-probA
- tempA <- sapply(probA, function(i) rbinom(grpA.size, ploidy, i) )
- if(grpA.size==1) {tempA <- matrix(tempA,nrow=1)}
- tempB <- sapply(probB, function(i) rbinom(grpB.size, ploidy, i) )
- if(grpB.size==1) {tempB <- matrix(tempB,nrow=1)}
- return(new("genlight", rbind(tempA,tempB), ploidy=ploidy))
- }
- }
-
- ## NON-STRUCTURED DATA ##
- ## generate data
- if(n.snp.nonstruc <= block.size){ # no need to use blocks
- res.ns <- f1(n.ind, n.snp.nonstruc)
- } else { # proceed by blocks of SNPs
- res.ns <- f1(n.ind, block.size)
- snp.to.go <- n.snp.nonstruc - nLoc(res.ns) # number of SNPs left to simulate
- while(snp.to.go > 0){
- nb.snp.tmp <- min(snp.to.go, block.size)
- res.ns <- cbind(res.ns, f1(n.ind, nb.snp.tmp))
- snp.to.go <- n.snp.nonstruc - nLoc(res.ns) # number of SNPs left to simulate
- }
- } # end non-structured data
-
-
- ## STRUCTURED DATA ##
- if(n.snp.struc > 0){
- if(n.snp.struc <= block.size){ # no need to use blocks
- res.s <- f2(n.ind, n.snp.struc)
- } else { # proceed by blocks of SNPs
- res.s <- f2(n.ind, block.size)
- snp.to.go <- n.snp.struc - nLoc(res.s) # number of SNPs left to simulate
- while(snp.to.go > 0){
- nb.snp <- min(snp.to.go, block.size)
- res.s <- cbind(res.s, f2(n.ind, nb.snp))
- snp.to.go <- n.snp.struc - nLoc(res.s) # number of SNPs left to simulate
- }
-
- }
- } # end structured data
-
- ## RETURN RESULTS ##
- res <- res.ns
- if(n.snp.struc > 0){
- res <- cbind(res,res.s)
- }
- if(!is.null(grp.size)){
- pop(res) <- rep(c("A","B"), c(grpA.size, grpB.size))
- }
- indNames(res) <- paste("ind", 1:n.ind)
-
- return(res)
-
-} # end glSim
-
-
-
-
-
-
###################
### TESTING
###################
Added: pkg/R/glSim.R
===================================================================
--- pkg/R/glSim.R (rev 0)
+++ pkg/R/glSim.R 2011-03-10 16:42:42 UTC (rev 849)
@@ -0,0 +1,96 @@
+
+##########
+## glSim
+##########
+glSim <- function(n.ind, n.snp.nonstruc, n.snp.struc=0, grp.size=round(n.ind/2), ploidy=1, alpha=0,
+ block.size=NULL){
+
+ ## BASIC CHECKS ##
+ if( any(c(n.ind, n.snp.nonstruc+n.snp.struc) <1)) stop("null numbers of individuals and/or SNPs requested")
+ ## alpha parameter
+ if(alpha>0.5){
+ alpha <- 0.5
+ warning("alpha cannot exceed 0.5 - changing the value to 0.5 (total forced asymmetry)")
+ }
+ if(alpha<0){
+ alpha <- 0
+ warning("alpha cannot be lower than 0 - changing the value to 0 (no forced asymmetry)")
+ }
+
+ ## handle block size
+ if(is.null(block.size)){
+ block.size <- n.snp.nonstruc + n.snp.struc
+ }
+
+ ## handle group sizes
+ grpA.size <- grp.size
+ if(grpA.size >= n.ind) stop("grpA.size is >= n.ind")
+ grpB.size <- n.ind - grpA.size
+
+
+ ## AUXIL FUNCTIONS ##
+ ## draw p snp for i indiv and convert into a genlight - no structure
+ f1 <- function(n,p){
+ temp <- sapply(1:p, function(i) rbinom(n, ploidy, runif(1)))
+ if(n==1) {temp <- matrix(temp,nrow=1)}
+ return(new("genlight", temp, ploidy=ploidy))
+ }
+
+ ## draw p snp for i indiv and convert into a genlight - differences between 2 groups
+ if(n.snp.struc > 0){
+ f2 <- function(n,p){
+ probA <- runif(p, min=0, max=0.5-alpha)
+ probB <- 1-probA
+ tempA <- sapply(probA, function(i) rbinom(grpA.size, ploidy, i) )
+ if(grpA.size==1) {tempA <- matrix(tempA,nrow=1)}
+ tempB <- sapply(probB, function(i) rbinom(grpB.size, ploidy, i) )
+ if(grpB.size==1) {tempB <- matrix(tempB,nrow=1)}
+ return(new("genlight", rbind(tempA,tempB), ploidy=ploidy))
+ }
+ }
+
+ ## NON-STRUCTURED DATA ##
+ ## generate data
+ if(n.snp.nonstruc <= block.size){ # no need to use blocks
+ res.ns <- f1(n.ind, n.snp.nonstruc)
+ } else { # proceed by blocks of SNPs
+ res.ns <- f1(n.ind, block.size)
+ snp.to.go <- n.snp.nonstruc - nLoc(res.ns) # number of SNPs left to simulate
+ while(snp.to.go > 0){
+ nb.snp.tmp <- min(snp.to.go, block.size)
+ res.ns <- cbind(res.ns, f1(n.ind, nb.snp.tmp))
+ snp.to.go <- n.snp.nonstruc - nLoc(res.ns) # number of SNPs left to simulate
+ }
+ } # end non-structured data
+
+
+ ## STRUCTURED DATA ##
+ if(n.snp.struc > 0){
+ if(n.snp.struc <= block.size){ # no need to use blocks
+ res.s <- f2(n.ind, n.snp.struc)
+ } else { # proceed by blocks of SNPs
+ res.s <- f2(n.ind, block.size)
+ snp.to.go <- n.snp.struc - nLoc(res.s) # number of SNPs left to simulate
+ while(snp.to.go > 0){
+ nb.snp <- min(snp.to.go, block.size)
+ res.s <- cbind(res.s, f2(n.ind, nb.snp))
+ snp.to.go <- n.snp.struc - nLoc(res.s) # number of SNPs left to simulate
+ }
+
+ }
+ } # end structured data
+
+ ## RETURN RESULTS ##
+ res <- res.ns
+ if(n.snp.struc > 0){
+ res <- cbind(res,res.s)
+ }
+ if(!is.null(grp.size)){
+ pop(res) <- rep(c("A","B"), c(grpA.size, grpB.size))
+ }
+ indNames(res) <- paste("ind", 1:n.ind)
+
+ return(res)
+
+} # end glSim
+
More information about the adegenet-commits
mailing list