diff --git a/.Rbuildignore b/.Rbuildignore new file mode 100644 index 0000000..91114bf --- /dev/null +++ b/.Rbuildignore @@ -0,0 +1,2 @@ +^.*\.Rproj$ +^\.Rproj\.user$ diff --git a/.gitignore b/.gitignore index c066a39..8e3a7fd 100644 --- a/.gitignore +++ b/.gitignore @@ -1 +1,2 @@ inst/doc +.Rproj.user diff --git a/DESCRIPTION b/DESCRIPTION index 30638bb..4500fff 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -8,5 +8,9 @@ Maintainer: Diptavo Dutta Description: Kernel Regression based association tests for Multiple phenotypes. The functions aggregate variant-phenotype score statistic in a particular region and computes corresponding p-values efficiently. Depends: SKAT, nlme, copula License: GPL (>=2) -Suggests: knitr,rmarkdown,R.rsp +Suggests: + knitr, + rmarkdown, + R.rsp, + testthat VignetteBuilder: knitr diff --git a/MultiSKAT.Rproj b/MultiSKAT.Rproj new file mode 100644 index 0000000..21a4da0 --- /dev/null +++ b/MultiSKAT.Rproj @@ -0,0 +1,17 @@ +Version: 1.0 + +RestoreWorkspace: Default +SaveWorkspace: Default +AlwaysSaveHistory: Default + +EnableCodeIndexing: Yes +UseSpacesForTab: Yes +NumSpacesForTab: 2 +Encoding: UTF-8 + +RnwWeave: Sweave +LaTeX: pdfLaTeX + +BuildType: Package +PackageUseDevtools: Yes +PackageInstallArgs: --no-multiarch --with-keep.source diff --git a/R/Genotype.Kernels.R b/R/Genotype.Kernels.R index c7a1204..ab78061 100644 --- a/R/Genotype.Kernels.R +++ b/R/Genotype.Kernels.R @@ -11,9 +11,9 @@ function(Z,obj.res, kernel = "linear.weighted", Is.Common=FALSE, weights.beta=c( if(m < 3 & verbose){ msg <-sprintf("Rare variant test with < 3 variants is not advisable") warning(msg,call. = FALSE)} - n.rare <- sum(MAF(Z) < 0.01) + n.rare <- sum(MAF(Z) < 0.01, na.rm = TRUE) if(verbose){print(paste("The region has ",n.rare," rare variants",sep = ""))} - mc <- sum(MAC(Z)) + mc <- sum(MAC(Z), na.rm = TRUE) if(mc < 5 & verbose){ msg <-sprintf("Rare variant test with total MAC < 5 is not advisable") warning(msg,call. = FALSE)} diff --git a/R/Helper.R b/R/Helper.R index 455272d..08e82e2 100644 --- a/R/Helper.R +++ b/R/Helper.R @@ -1,28 +1,28 @@ -mat.sqrt <- function(A) -{ - ei<-eigen(A) - d<-ei$values - d<-(d+abs(d))/2 - d2<-sqrt(d) - ans<-ei$vectors %*% diag(d2) %*% t(ei$vectors) - return(ans) -} - -MAF <- function(G){ - mf <- array() - for(i in 1:ncol(G)){ - mf[i] <- sqrt(sum(G[,i] == 0)/nrow(G)); - if(mf[i] > 0.5) - mf[i] = 1- mf[i]; - } - return(mf) -} - -MAC <- function(G){ - mc <- array() - for(i in 1:ncol(G)){ - mc[i] <- sum(G[,i] == 1) + 2*min(sum(G[,i] == 2),sum(G[,i] == 0)); - } - return(mc); -} - +mat.sqrt <- function(A) +{ + ei<-eigen(A) + d<-ei$values + d<-(d+abs(d))/2 + d2<-sqrt(d) + ans<-ei$vectors %*% diag(d2) %*% t(ei$vectors) + return(ans) +} + +MAF <- function(G){ + mf <- array() + for(i in 1:ncol(G)){ + mf[i] <- sqrt(sum(G[,i] == 0, na.rm = TRUE)/nrow(G)); + if(mf[i] > 0.5) + mf[i] = 1- mf[i]; + } + return(mf) +} + +MAC <- function(G){ + mc <- array() + for(i in 1:ncol(G)){ + mc[i] <- sum(G[,i] == 1, na.rm = TRUE) + 2*min(sum(G[,i] == 2, na.rm = TRUE), sum(G[,i] == 0, na.rm = TRUE), na.rm = TRUE); + } + return(mc); +} + diff --git a/tests/testthat.R b/tests/testthat.R new file mode 100644 index 0000000..abdc476 --- /dev/null +++ b/tests/testthat.R @@ -0,0 +1,4 @@ +library(testthat) +library(MultiSKAT) + +test_check("MultiSKAT") diff --git a/tests/testthat/test-missing-genotype-data.R b/tests/testthat/test-missing-genotype-data.R new file mode 100644 index 0000000..20aa1c4 --- /dev/null +++ b/tests/testthat/test-missing-genotype-data.R @@ -0,0 +1,18 @@ +context("test-missing-genotype-data") + + + +test_that("multiplication works", { + + data(MultiSKAT.example.data) + attach(MultiSKAT.example.data) + + ## Introduce NAs in Genotypes + Genotypes[1,1] <- NA + + expect_warning(MultiSKAT(MultiSKAT_NULL(Phenotypes, Cov), Genotypes, Sigma_p = cov(Phenotypes), verbose = FALSE), + "The missing genotype rate is 0.000004. Imputation is applied.") + + detach(MultiSKAT.example.data) + +})