Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Enable imputed dosage for admixmap() #119

Open
wants to merge 2 commits into
base: devel
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions .Rbuildignore
Original file line number Diff line number Diff line change
@@ -1,3 +1,5 @@
^renv$
^renv\.lock$
^\.github$
^.*\.Rproj$
^\.Rproj\.user$
1 change: 1 addition & 0 deletions GENESIS
Submodule GENESIS added at dd2a91
13 changes: 11 additions & 2 deletions R/admixMap.R
Original file line number Diff line number Diff line change
@@ -1,6 +1,7 @@
admixMap <- function(admixDataList,
null.model,
male.diploid=TRUE, genome.build=c("hg19", "hg38"),
male.diploid=TRUE, genome.build=c("hg19", "hg38"),imputed=FALSE,
dosage.field="DS",
BPPARAM=bpparam(), verbose=TRUE){

# if admixDataList is one file, convert to a list
Expand All @@ -15,7 +16,7 @@ admixMap <- function(admixDataList,
if(is.null(names(admixDataList))){
names(admixDataList) <- paste("Anc",1:v,sep="")
}

# get sample index
if (is(admixDataList[[1]], "GenotypeIterator")) {
sample.index <- lapply(admixDataList, .sampleIndexNullModel, null.model)
Expand All @@ -30,6 +31,10 @@ admixMap <- function(admixDataList,
}
n.samp <- length(sample.index)

# if admixDataList contain GenotypeIterator, imputed=TRUE will be ignored
if (is(admixDataList[[1]], "GenotypeIterator") && imputed=TRUE)
print("admixDataList contain GenotypeIterator, imputed=TRUE is ignored")

# get sex for calculating allele freq
sex <- validateSex(admixDataList[[1]])[sample.index]

Expand Down Expand Up @@ -81,7 +86,11 @@ admixMap <- function(admixDataList,
if (is(admixDataList[[1]], "GenotypeIterator")) {
local[,,i] <- getGenotypeSelection(admixDataList[[i]], scan=sample.index, order="selection", transpose=TRUE, use.names=FALSE, drop=FALSE)
} else {
if(imputed=FALSE){
local[,,i] <- refDosage(admixDataList[[i]], use.names=FALSE)[sample.index,,drop=FALSE]
}else
local[,,i] <- imputedDosage(admixDataList[[i]], use.names=FALSE,dosage.field="DS")[sample.index,,drop=FALSE]

}
}
if (any(is.na(local))) warning("missing values in local ancestry will produce NA output for this block")
Expand Down
5 changes: 4 additions & 1 deletion man/admixMap.Rd
Original file line number Diff line number Diff line change
Expand Up @@ -7,7 +7,8 @@ Run admixture analyses
}
\usage{
admixMap(admixDataList, null.model, male.diploid=TRUE,
genome.build=c("hg19", "hg38"),
genome.build=c("hg19", "hg38"),imputed=FALSE,
dosage.field="DS",
BPPARAM=bpparam(), verbose=TRUE)
}

Expand All @@ -16,6 +17,8 @@ admixMap(admixDataList, null.model, male.diploid=TRUE,
\item{null.model}{A null model object returned by \code{\link{fitNullModel}}.}
\item{male.diploid}{Logical for whether males on sex chromosomes are coded as diploid. Default is `male.diploid=TRUE`, meaning sex chromosome genotypes for males have values 0/2. If the input object codes males as 0/1 on sex chromosomes, set `male.diploid=FALSE`.}
\item{genome.build}{A character sting indicating genome build; used to identify pseudoautosomal regions on the X and Y chromosomes. These regions are not treated as sex chromosomes when calculating allele frequencies.}
\item{imputed}{Logical indicator of whether to read dosages from the field specified in `dosage.field`,containing imputed dosages instead of counting the number of ref alleles. Only relevant when admixDataList contains \code{\link{SeqVarIterator}}, i.e. ignored when admixDataList contains \code{\link{GenotypeIterator}} }
\item{dosage.field}{The name of the dosage field in the GDS object (will be prepended with "annotation/format"). Default is dosage.field='DS'.}
\item{BPPARAM}{A \code{\link{BiocParallelParam}} object to process blocks of variants in parallel. If not provided, the default back-end returned by \code{\link{bpparam}} will be used.}
\item{verbose}{Logical indicator of whether updates from the function should be printed to the console; the default is TRUE.}
}
Expand Down