@@ -9,7 +9,7 @@ Description: This package implements specialized algorithms that enable

     following publication: Belleau, P et al. Genetic Ancestry Inference from 

     Cancer-Derived Molecular Data across Genomic and Transcriptomic 

     Platforms. Cancer Res 1 January 2023; 83 (1): 49–58.

-Version: 1.3.1

+Version: 1.3.2

 Authors@R: c(person("Pascal", "Belleau", email="[email protected]",

     role=c("cre", "aut"), comment = c(ORCID = "0000-0002-0802-1071")),

     person("Astrid", "Deschênes", email="[email protected]",

@@ -1,6 +1,7 @@

 # Generated by roxygen2: do not edit by hand

 export(add1KG2SampleGDS)

+export(addBlockFromDetFile)

 export(addGeneBlockGDSRefAnnot)

 export(addGeneBlockRefAnnot)

 export(addRef2GDS1KG)

@@ -1335,3 +1335,107 @@ appendGDSSampleOnly <- function(gds, listSamples) {

     return(0L)

 }

+

+#' @title Append information associated to ld blocks, as indexes, into the

+#' Population Reference SNV Annotation GDS file

+#'

+#' @description The function appends the information about the ld blocks into

+#' the Population Reference SNV Annotation GDS file. The information is

+#' extracted from the parameter listBlock.

+#'

+#' @param gds an object of class \link[gdsfmt]{gds.class}

+#' (GDS file), an opened Reference Annotation GDS file.

+#'

+#' @param listBlock a \code{array} of integer

+#' representing the linkage disequilibrium block for

+#' each SNV in the in the same order than the variant

+#' in Population reference dataset.

+#'

+#' @param blockName a \code{character} string representing the id of the block.

+#' The blockName should not exist in \'gdsRefAnnotFile\'.

+#'

+#' @param blockDesc a \code{character} string representing the description of

+#' the block.

+#'

+#' @return The integer \code{0L} when successful.

+#'

+#' @examples

+#'

+#' ## Required library for GDS

+#' library(gdsfmt)

+#' ## Path to the demo pedigree file is located in this package

+#' dataDir <- system.file("extdata", package="RAIDS")

+#'

+#  ## Temporary file

+#' fileAnnotGDS <- file.path(tempdir(), "ex1_good_small_1KG_Ann_GDS.gds")

+#'

+#'

+#' file.copy(file.path(dataDir, "tests",

+#'     "ex1_NoBlockGene.1KG_Annot_GDS.gds"), fileAnnotGDS)

+#'

+#'

+#' fileReferenceGDS  <- file.path(dataDir, "tests",

+#'     "ex1_good_small_1KG.gds")

+#'  \donttest{

+#'     gdsRef <- openfn.gds(fileReferenceGDS)

+#'     listBlock <- read.gdsn(index.gdsn(gdsRef, "snp.position"))

+#'     listBlock <- rep(-1, length(listBlock))

+#'     closefn.gds(gdsRef)

+#'     gdsAnnot1KG <- openfn.gds(fileAnnotGDS, readonly=FALSE)

+#'     ## Append information associated to blocks

+#'     RAIDS:::addGDS1KGLDBlock(gds=gdsAnnot1KG,

+#'         listBlock=listBlock,

+#'         blockName="blockEmpty",

+#'         blockDesc="Example")

+#'.    closefn.gds(gdsAnnot1KG)

+#'

+#'     gdsAnnot1KG <- openfn.gds(fileAnnotGDS)

+#'     print(gdsAnnot1KG)

+#'

+#'     closefn.gds(gdsAnnot1KG)

+#' }

+#'

+#' ## Remove temporary file

+#' unlink(fileAnnotGDS, force=TRUE)

+#'

+#' @author Pascal Belleau, Astrid Deschênes and Alexander Krasnitz

+#' @importFrom gdsfmt add.gdsn index.gdsn ls.gdsn compression.gdsn

+#' @importFrom gdsfmt append.gdsn sync.gds

+#' @encoding UTF-8

+#' @keywords internal

+addGDS1KGLDBlock <- function(gds, listBlock, blockName, blockDesc) {

+

+    blockAnnot <- data.frame(block.id=blockName,

+                             block.desc=blockDesc,

+                             stringsAsFactors=FALSE)

+

+    if(! ("block.annot" %in% ls.gdsn(gds))) {

+        varBlockAnnot <- add.gdsn(gds, "block.annot", blockAnnot)

+    }else {

+        curAnnot <- index.gdsn(gds, "block.annot/block.id")

+        append.gdsn(curAnnot,blockAnnot$block.id)

+        curAnnot <- index.gdsn(gds, "block.annot/block.desc")

+        append.gdsn(curAnnot, blockAnnot$block.desc)

+    }

+

+    varBlock <- NULL

+    if(! ("block" %in% ls.gdsn(gds))){

+        varBlock <- add.gdsn(gds, "block",

+                             valdim=c(length(listBlock), 1),

+                             listBlock, storage="int32",

+                             compress = "LZ4_RA")

+        readmode.gdsn(varBlock)

+

+    }else {

+        if(is.null(varBlock)) {

+            varBlock <- index.gdsn(gds, "block")

+            varBlock <- compression.gdsn(varBlock, "")

+        }

+        append.gdsn(varBlock, listBlock)

+        varBlock <- compression.gdsn(varBlock, "LZ4_RA")

+    }

+

+    sync.gds(gds)

+

+    return(0L)

+}

@@ -1020,6 +1020,183 @@ getRefSuperPop <- function(fileReferenceGDS) {

     return(df)

 }

+#' @title Append information associated to ld blocks, as indexes, into the

+#' Population Reference SNV Annotation GDS file

+#'

+#' @description The function appends the information about the ld blocks into

+#' the Population Reference SNV Annotation GDS file. The information is

+#' extracted from the Population Reference GDS file and files \'.det\'.

+#'

+#' @param fileReferenceGDS a \code{character} string representing the file

+#' name of the Reference GDS file. The file must exist.

+#'

+#' @param gdsRefAnnotFile a \code{character} string representing the

+#' file name corresponding the Reference SNV

+#' Annotation GDS file. The function will

+#' open it in write mode and close it after. The file must exist.

+#'

+#' @param pathBlock a \code{character} string representing the directory

+#' where all the output file det from the plink block command are located.

+#' The directory must not include other file with the extension \'.det\'.

+#' The name of the \'.det\' must include the super-population between \'.\'

+#' and the chromosome in the form \'chrNumber.\' \( \'chr1.\'\).

+#'

+#' @param superPop a \code{character} string representing the super population.

+#'

+#' @param blockName a \code{character} string representing the id of the block.

+#' The blockName should not exist in \'gdsRefAnnotFile\'.

+#' Default: \code{"ldBlock"}.

+#'

+#' @param blockDesc a \code{character} string representing the description of

+#' the block.

+#' Default: \code{"Not Define"}

+#'

+#' @param verbose a \code{logical} indicating if message information should be

+#' printed. Default: \code{FALSE}.

+#'

+#' @return \code{OL} when the function is successful.

+#'

+#' @details

+#'

+#' More information about GDS file format can be found at the Bioconductor

+#' gdsfmt website:

+#' https://bioconductor.org/packages/gdsfmt/

+#'

+#' @examples

+#'

+#' ## Path to the demo pedigree file is located in this package

+#' dataDir <- system.file("extdata", package="RAIDS")

+#'

+#  ## Temporary file

+#' fileAnnotGDS <- file.path(tempdir(), "ex1_good_small_1KG_Ann_GDS.gds")

+#'

+#' ## Demo of of output file det from the plink block

+#' ## command for chromosome 1

+#' fileLdBlock <- file.path(dirname(fileAnnotGDS), "block.sp.EUR.Ex.chr1.blocks.det")

+#'

+#'

+#' file.copy(file.path(dataDir, "tests",

+#'     "ex1_NoBlockGene.1KG_Annot_GDS.gds"), fileAnnotGDS)

+#' file.copy(file.path(dataDir, "block.sp.EUR.Ex.chr1.blocks.det"),

+#'     fileLdBlock)

+#'

+#'

+#'

+#' ## GDS Reference file

+#' fileReferenceGDS  <- file.path(dataDir, "tests",

+#'     "ex1_good_small_1KG.gds")

+#'

+#'  \donttest{

+#'

+#'

+#'     ## Append information associated to blocks

+#'     addBlockFromDetFile(fileReferenceGDS=fileReferenceGDS,

+#'         gdsRefAnnotFile=fileAnnotGDS,

+#'         pathBlock=dirname(fileAnnotGDS),

+#'         superPop="EUR")

+#'

+#'     gdsAnnot1KG <- openfn.gds(fileAnnotGDS)

+#'     print(gdsAnnot1KG)

+#'

+#'     closefn.gds(gdsAnnot1KG)

+#' }

+#'

+#' ## Remove temporary file

+#' unlink(fileAnnotGDS, force=TRUE)

+#' unlink(fileLdBlock, force=TRUE)

+#'

+#' @author Pascal Belleau, Astrid Deschênes and Alexander Krasnitz

+#'

+#' @importFrom gdsfmt openfn.gds closefn.gds read.gdsn index.gdsn ls.gdsn

+#' @importFrom SNPRelate snpgdsOpen

+#' @encoding UTF-8

+#' @export

+addBlockFromDetFile <- function(fileReferenceGDS, gdsRefAnnotFile, pathBlock,

+                                superPop, blockName="ldBlock",

+                                blockDesc="Not Define", verbose=FALSE) {

+    if (!(is.character(fileReferenceGDS) && (file.exists(fileReferenceGDS)))) {

+        stop("The \'fileReferenceGDS\' must be a character string ",

+             "representing the Reference GDS file. The file must exist.")

+    }

+    if(!(is.character(blockName))){

+        stop("The \'blockName\' must be a character string ",

+             "representing the name of the block.")

+    }

+

+    if(blockName == "ldBlock"){

+        blockName <- paste0(blockName, ".", superPop)

+    }

+

+    gdsRefAnnot <- openfn.gds(gdsRefAnnotFile)

+

+    if(("block.annot" %in% ls.gdsn(gdsRefAnnot))) {

+        listAnno <- read.gdsn(index.gdsn(gdsRefAnnot, "block.annot"))

+        if(length(which(gdsRefAnnot$block.id  == blockName)) > 0){

+            stop("The \'blockName\' already exist in \'gdsRefAnnotFile\'.")

+        }

+    }

+    closefn.gds(gdsRefAnnot)

+

+    gdsReference <- snpgdsOpen(filename=fileReferenceGDS)

+

+

+

+    ## The verbose must be a logical

+    validateLogical(verbose, "verbose")

+

+    ## Extract the SNP chromosomes and positions

+    snpChromosome <- read.gdsn(index.gdsn(gdsReference, "snp.chromosome"))

+    #snpPosition <- read.gdsn(index.gdsn(gdsReference, "snp.position"))

+    closefn.gds(gdsReference)

+

+    listFileBlock <- dir(pathBlock, ".det")

+    listFileBlock <- listFileBlock[grep(paste0("\\.", superPop, "\\."), listFileBlock)]

+

+    listChr <- unique(snpChromosome)

+

+    #listChr <- listChr[order(listChr)]

+    #listChr <- seq_len(22)

+    listBlock <- list()

+

+    for(chr in seq_len(length(listChr))) {

+        if(verbose) { message("chr", listChr[chr], " ",Sys.time()) }

+        listChrCur <- listFileBlock[grep(paste0("chr",listChr[chr],"\\."), listFileBlock)]

+        if(length(listChrCur) == 1){

+            tmp <- processBlockChr(fileReferenceGDS, file.path(pathBlock, listChrCur))

+            listBlock[[chr]] <- tmp$block.snp

+            if(chr > 1) {

+                vMax <- max(listBlock[[chr-1]])

+                vMin <- min(listBlock[[chr-1]])

+                listBlock[[chr]][listBlock[[chr]] > 0] <-

+                    listBlock[[chr]][listBlock[[chr]] > 0] + vMax

+                if(vMin < 0) {

+                    listBlock[[chr]][listBlock[[chr]] < 0] <-

+                        listBlock[[chr]][listBlock[[chr]] < 0] + vMin

+                }

+            }

+        }else{

+

+            listBlock[[chr]] <- rep(-1, length(which(snpChromosome == listChr[chr])))

+            vMin <- 0

+            if(chr > 1){

+                vMin <- min(listBlock[[chr-1]])

+            }

+            if(vMin < 0){

+                listBlock[[chr]] <- listBlock[[chr]] + vMin

+            }

+        }

+

+    }

+    listBlock <- do.call(c, listBlock)

+

+    gdsRefAnnot <- openfn.gds(gdsRefAnnotFile, readonly = FALSE)

+

+    ## Save the information into the GDS file

+    addGDS1KGLDBlock(gdsRefAnnot, listBlock, blockName, blockDesc)

+    closefn.gds(gdsRefAnnot)

+    ## Success

+    return(0L)

+}

 #' @title Append information associated to blocks, as indexes, into the

 #' Population Reference SNV Annotation GDS file

@@ -2351,7 +2351,9 @@ inferAncestry <- function(profileFile, pathProfileGDS,

     genoSource <- arg_match(genoSource)

-

+    if(genoSource == "bam"){

+        stop("The bam is not release yet look to get a \'Devel\' version or contact us")

+    }

     profileName <- gsub("\\.gz$", "", profileBaseName, ignore.case = TRUE)

     for(extCur in c( "\\.vcf$", "\\.txt$", "\\.bam", "\\.tsv", "\\.csv")){

         profileName <- gsub(extCur, "", profileName, ignore.case = TRUE)

@@ -2776,6 +2778,9 @@ inferAncestryGeneAware <- function(profileFile, pathProfileGDS,

     genoSource <- arg_match(genoSource)

+    if(genoSource == "bam"){

+        stop("The bam is not release yet look to get a \'Devel\' version or contact us")

+    }

     profileName <- gsub("\\.gz$", "", profileBaseName, ignore.case = TRUE)

     for(extCur in c( "\\.vcf$", "\\.txt$", "\\.bam", "\\.tsv", "\\.csv")){

@@ -458,7 +458,7 @@ readSNVVCF <- function(fileName,

 #' \item{\code{chr}}{ a \code{integer} representing a the chromosome from

 #' fileBlock.

 #' }

-#' \item{\code{block.snp}}{ the a \code{array} of integer

+#' \item{\code{block.snp}}{ a \code{array} of integer

 #' representing the linkage disequilibrium block for

 #' each SNV in the in the same order than the variant

 #' in Population reference dataset.

@@ -506,43 +506,46 @@ processBlockChr <- function(fileReferenceGDS, fileBlock) {

     listChr <- as.integer(gsub("chr", "", listChr))

     listSNVChr <- read.gdsn(index.gdsn(gdsReference, "snp.chromosome"))

     listSNVChr <- which(listSNVChr == listChr)

-    snp.keep <- read.gdsn(index.gdsn(gdsReference, "snp.position"))[listSNVChr]

+    snpKeep <- read.gdsn(index.gdsn(gdsReference, "snp.position"))[listSNVChr]

     closefn.gds(gdsReference)

-    z <- cbind(c(blockChr$BP1, snp.keep, blockChr$BP2+1),

+    z <- cbind(c(blockChr$BP1, snpKeep, blockChr$BP2+1),

                c(seq_len(nrow(blockChr)),

-                 rep(0, length(snp.keep)), -1*seq_len(nrow(blockChr))))

+                 rep(0, length(snpKeep)), -1*seq_len(nrow(blockChr))))

     z <- z[order(z[,1]),]

-    block.snp <- cumsum(z[,2])[z[,2] == 0]

+    blockSnp <- cumsum(z[,2])[z[,2] == 0]

     curStart <- 0

     activeBlock <- 0

     blockState <- 0

-    block.inter <- rep(0, length(which(block.snp == 0)))

+    blockInter <- rep(0, length(which(blockSnp == 0)))

     k <- 1

-    for(i in seq_len(length(block.snp))){

-        if(block.snp[i] == 0){

+    for(i in seq_len(length(blockSnp))){

+        if(blockSnp[i] == 0){

             if(activeBlock == 1){

-                if(snp.keep[i] - curStart >= 10000) {

+                if(snpKeep[i] - curStart >= 10000) {

                     blockState <- blockState - 1

-                    curStart <- snp.keep[i]

+                    curStart <- snpKeep[i]

                 }

             } else{

                 blockState <- blockState - 1

-                curStart <- snp.keep[i]

-                curStart <- snp.keep[i]

+                curStart <- snpKeep[i]

                 activeBlock <- 1

             }

-            block.inter[k] <- blockState

+            if(blockState == 0){

+                blockState <- -1

+            }

+            blockInter[k] <- blockState

             k <- k + 1

         }else{

             activeBlock <- 0

         }

     }

-    block.snp[block.snp == 0] <- block.inter

+

+    blockSnp[blockSnp == 0] <- blockInter

     res <- list(chr=listChr,

-                block.snp=block.snp)

+                block.snp=blockSnp)

     return(res)

 }

@@ -1,3 +1,10 @@

+CHANGES IN VERSION 1.3.2

+------------------------

+

+SIGNIFICANT USER-VISIBLE CHANGES

+

+    o New functions inferAncestry(), inferAncestryGeneAware() and getRefSuperPop() to simplify ancestry inference

+

 CHANGES IN VERSION 1.3.1

 ------------------------

new file mode 100644

@@ -0,0 +1,105 @@

+% Generated by roxygen2: do not edit by hand

+% Please edit documentation in R/process1KG.R

+\encoding{UTF-8}

+\name{addBlockFromDetFile}

+\alias{addBlockFromDetFile}

+\title{Append information associated to ld blocks, as indexes, into the

+Population Reference SNV Annotation GDS file}

+\usage{

+addBlockFromDetFile(

+  fileReferenceGDS,

+  gdsRefAnnotFile,

+  pathBlock,

+  superPop,

+  blockName = "ldBlock",

+  blockDesc = "Not Define",

+  verbose = FALSE

+)

+}

+\arguments{

+\item{fileReferenceGDS}{a \code{character} string representing the file

+name of the Reference GDS file. The file must exist.}

+

+\item{gdsRefAnnotFile}{a \code{character} string representing the

+file name corresponding the Reference SNV

+Annotation GDS file. The function will

+open it in write mode and close it after. The file must exist.}

+

+\item{pathBlock}{a \code{character} string representing the directory

+where all the output file det from the plink block command are located.

+The directory must not include other file with the extension \'.det\'.

+The name of the \'.det\' must include the super-population between \'.\'

+and the chromosome in the form \'chrNumber.\' \( \'chr1.\'\).}

+

+\item{superPop}{a \code{character} string representing the super population.}

+

+\item{blockName}{a \code{character} string representing the id of the block.

+The blockName should not exist in \'gdsRefAnnotFile\'.

+Default: \code{"ldBlock"}.}

+

+\item{blockDesc}{a \code{character} string representing the description of

+the block.

+Default: \code{"Not Define"}}

+

+\item{verbose}{a \code{logical} indicating if message information should be

+printed. Default: \code{FALSE}.}

+}

+\value{

+\code{OL} when the function is successful.

+}

+\description{

+The function appends the information about the ld blocks into

+the Population Reference SNV Annotation GDS file. The information is

+extracted from the Population Reference GDS file and files \'.det\'.

+}

+\details{

+More information about GDS file format can be found at the Bioconductor

+gdsfmt website:

+https://bioconductor.org/packages/gdsfmt/

+}

+\examples{

+

+## Path to the demo pedigree file is located in this package

+dataDir <- system.file("extdata", package="RAIDS")

+

+fileAnnotGDS <- file.path(tempdir(), "ex1_good_small_1KG_Ann_GDS.gds")

+

+## Demo of of output file det from the plink block

+## command for chromosome 1

+fileLdBlock <- file.path(dirname(fileAnnotGDS), "block.sp.EUR.Ex.chr1.blocks.det")

+

+

+file.copy(file.path(dataDir, "tests",

+    "ex1_NoBlockGene.1KG_Annot_GDS.gds"), fileAnnotGDS)

+file.copy(file.path(dataDir, "block.sp.EUR.Ex.chr1.blocks.det"),

+    fileLdBlock)

+

+

+

+## GDS Reference file

+fileReferenceGDS  <- file.path(dataDir, "tests",

+    "ex1_good_small_1KG.gds")

+

+ \donttest{

+

+

+    ## Append information associated to blocks

+    addBlockFromDetFile(fileReferenceGDS=fileReferenceGDS,

+        gdsRefAnnotFile=fileAnnotGDS,

+        pathBlock=dirname(fileAnnotGDS),

+        superPop="EUR")

+

+    gdsAnnot1KG <- openfn.gds(fileAnnotGDS)

+    print(gdsAnnot1KG)

+

+    closefn.gds(gdsAnnot1KG)

+}

+

+## Remove temporary file

+unlink(fileAnnotGDS, force=TRUE)

+unlink(fileLdBlock, force=TRUE)

+

+}

+\author{

+Pascal Belleau, Astrid Deschênes and Alexander Krasnitz

+}

new file mode 100644

@@ -0,0 +1,76 @@

+% Generated by roxygen2: do not edit by hand

+% Please edit documentation in R/gdsWrapper_internal.R

+\encoding{UTF-8}

+\name{addGDS1KGLDBlock}

+\alias{addGDS1KGLDBlock}

+\title{Append information associated to ld blocks, as indexes, into the

+Population Reference SNV Annotation GDS file}

+\usage{

+addGDS1KGLDBlock(gds, listBlock, blockName, blockDesc)

+}

+\arguments{

+\item{gds}{an object of class \link[gdsfmt]{gds.class}

+(GDS file), an opened Reference Annotation GDS file.}

+

+\item{listBlock}{a \code{array} of integer

+representing the linkage disequilibrium block for

+each SNV in the in the same order than the variant

+in Population reference dataset.}

+

+\item{blockName}{a \code{character} string representing the id of the block.

+The blockName should not exist in \'gdsRefAnnotFile\'.}

+

+\item{blockDesc}{a \code{character} string representing the description of

+the block.}

+}

+\value{

+The integer \code{0L} when successful.

+}

+\description{

+The function appends the information about the ld blocks into

+the Population Reference SNV Annotation GDS file. The information is

+extracted from the parameter listBlock.

+}

+\examples{

+

+## Required library for GDS

+library(gdsfmt)

+## Path to the demo pedigree file is located in this package

+dataDir <- system.file("extdata", package="RAIDS")

+

+fileAnnotGDS <- file.path(tempdir(), "ex1_good_small_1KG_Ann_GDS.gds")

+

+

+file.copy(file.path(dataDir, "tests",

+    "ex1_NoBlockGene.1KG_Annot_GDS.gds"), fileAnnotGDS)

+

+

+fileReferenceGDS  <- file.path(dataDir, "tests",

+    "ex1_good_small_1KG.gds")

+ \donttest{

+    gdsRef <- openfn.gds(fileReferenceGDS)

+    listBlock <- read.gdsn(index.gdsn(gdsRef, "snp.position"))

+    listBlock <- rep(-1, length(listBlock))

+    closefn.gds(gdsRef)

+    gdsAnnot1KG <- openfn.gds(fileAnnotGDS, readonly=FALSE)

+    ## Append information associated to blocks

+    RAIDS:::addGDS1KGLDBlock(gds=gdsAnnot1KG,

+        listBlock=listBlock,

+        blockName="blockEmpty",

+        blockDesc="Example")

+.    closefn.gds(gdsAnnot1KG)

+

+    gdsAnnot1KG <- openfn.gds(fileAnnotGDS)

+    print(gdsAnnot1KG)

+

+    closefn.gds(gdsAnnot1KG)

+}

+

+## Remove temporary file

+unlink(fileAnnotGDS, force=TRUE)

+

+}

+\author{

+Pascal Belleau, Astrid Deschênes and Alexander Krasnitz

+}

+\keyword{internal}

@@ -22,7 +22,7 @@ a \code{list} containing 2 entries:

 \item{\code{chr}}{ a \code{integer} representing a the chromosome from

 fileBlock.

 }

-\item{\code{block.snp}}{ the a \code{array} of integer

+\item{\code{block.snp}}{ a \code{array} of integer

 representing the linkage disequilibrium block for

 each SNV in the in the same order than the variant

 in Population reference dataset.