Back to Mac ARM64 build report for BioC 3.17 |
|
This page was generated on 2023-10-20 09:37:58 -0400 (Fri, 20 Oct 2023).
Hostname | OS | Arch (*) | R version | Installed pkgs |
---|---|---|---|---|
kjohnson2 | macOS 12.6.1 Monterey | arm64 | 4.3.1 (2023-06-16) -- "Beagle Scouts" | 4347 |
Click on any hostname to see more info about the system (e.g. compilers) (*) as reported by 'uname -p', except on Windows and Mac OS X |
Package 109/2230 | Hostname | OS / Arch | INSTALL | BUILD | CHECK | BUILD BIN | ||||||||
atSNP 1.16.0 (landing page) Sunyoung Shin
| kjohnson2 | macOS 12.6.1 Monterey / arm64 | OK | OK | OK | OK | ||||||||
To the developers/maintainers of the atSNP package: - Use the following Renviron settings to reproduce errors and warnings. - If 'R CMD check' started to fail recently on the Linux builder(s) over a missing dependency, add the missing dependency to 'Suggests:' in your DESCRIPTION file. See Renviron.bioc for more information. |
Package: atSNP |
Version: 1.16.0 |
Command: /Library/Frameworks/R.framework/Resources/bin/R CMD check --install=check:atSNP.install-out.txt --library=/Library/Frameworks/R.framework/Resources/library --no-vignettes --timings atSNP_1.16.0.tar.gz |
StartedAt: 2023-10-17 03:25:17 -0400 (Tue, 17 Oct 2023) |
EndedAt: 2023-10-17 03:34:59 -0400 (Tue, 17 Oct 2023) |
EllapsedTime: 582.7 seconds |
RetCode: 0 |
Status: OK |
CheckDir: atSNP.Rcheck |
Warnings: 0 |
############################################################################## ############################################################################## ### ### Running command: ### ### /Library/Frameworks/R.framework/Resources/bin/R CMD check --install=check:atSNP.install-out.txt --library=/Library/Frameworks/R.framework/Resources/library --no-vignettes --timings atSNP_1.16.0.tar.gz ### ############################################################################## ############################################################################## * using log directory ‘/Users/biocbuild/bbs-3.17-bioc-mac-arm64/meat/atSNP.Rcheck’ * using R version 4.3.1 (2023-06-16) * using platform: aarch64-apple-darwin20 (64-bit) * R was compiled by Apple clang version 14.0.0 (clang-1400.0.29.202) GNU Fortran (GCC) 12.2.0 * running under: macOS Monterey 12.6.7 * using session charset: UTF-8 * using option ‘--no-vignettes’ * checking for file ‘atSNP/DESCRIPTION’ ... OK * checking extension type ... Package * this is package ‘atSNP’ version ‘1.16.0’ * checking package namespace information ... OK * checking package dependencies ... OK * checking if this is a source package ... OK * checking if there is a namespace ... OK * checking for hidden files and directories ... OK * checking for portable file names ... OK * checking for sufficient/correct file permissions ... OK * checking whether package ‘atSNP’ can be installed ... OK * used C++ compiler: ‘Apple clang version 14.0.0 (clang-1400.0.29.202)’ * used SDK: ‘MacOSX11.3.sdk’ * checking installed package size ... OK * checking package directory ... OK * checking ‘build’ directory ... OK * checking DESCRIPTION meta-information ... OK * checking top-level files ... OK * checking for left-over files ... OK * checking index information ... OK * checking package subdirectories ... OK * checking R files for non-ASCII characters ... OK * checking R files for syntax errors ... OK * checking whether the package can be loaded ... OK * checking whether the package can be loaded with stated dependencies ... OK * checking whether the package can be unloaded cleanly ... OK * checking whether the namespace can be loaded with stated dependencies ... OK * checking whether the namespace can be unloaded cleanly ... OK * checking dependencies in R code ... NOTE Namespaces in Imports field not imported from: ‘graphics’ ‘testthat’ All declared Imports should be used. * checking S3 generic/method consistency ... OK * checking replacement functions ... OK * checking foreign function calls ... OK * checking R code for possible problems ... NOTE ComputePValues: no visible binding for global variable ‘motif’ ComputePValues: no visible binding for global variable ‘snpid’ ComputePValues: no visible binding for global variable ‘snpbase’ ComputePValues: no visible binding for global variable ‘pval_ref’ ComputePValues: no visible binding for global variable ‘pval_snp’ ComputePValues: no visible binding for global variable ‘pval_cond_ref’ ComputePValues: no visible binding for global variable ‘pval_cond_snp’ ComputePValues: no visible binding for global variable ‘pval_diff’ ComputePValues: no visible binding for global variable ‘pval_rank’ LoadSNPData: no visible global function definition for ‘is’ LoadSNPData: no visible binding for global variable ‘IUPAC_CODE_MAP’ MatchSubsequence: no visible binding for global variable ‘motif’ MatchSubsequence: no visible binding for global variable ‘snpid’ MatchSubsequence: no visible binding for global variable ‘snpbase’ MatchSubsequence: no visible binding for global variable ‘len_seq’ MatchSubsequence: no visible binding for global variable ‘ref_seq’ MatchSubsequence : <anonymous>: no visible binding for global variable ‘motif’ checkMotifs: no visible global function definition for ‘is’ checkSNPids: no visible global function definition for ‘is’ dtMotifMatch: no visible binding for global variable ‘ref_seq’ dtMotifMatch: no visible binding for global variable ‘len_seq’ dtMotifMatch: no visible binding for global variable ‘snp_ref_start’ dtMotifMatch: no visible binding for global variable ‘ref_start’ dtMotifMatch: no visible binding for global variable ‘snp_start’ dtMotifMatch: no visible binding for global variable ‘snp_ref_end’ dtMotifMatch: no visible binding for global variable ‘ref_end’ dtMotifMatch: no visible binding for global variable ‘snp_end’ dtMotifMatch: no visible binding for global variable ‘snp_ref_length’ dtMotifMatch: no visible binding for global variable ‘ref_aug_match_seq_forward’ dtMotifMatch: no visible binding for global variable ‘ref_aug_match_seq_reverse’ dtMotifMatch: no visible binding for global variable ‘snp_aug_match_seq_forward’ dtMotifMatch: no visible binding for global variable ‘snp_seq’ dtMotifMatch: no visible binding for global variable ‘snp_aug_match_seq_reverse’ dtMotifMatch: no visible binding for global variable ‘ref_strand’ dtMotifMatch: no visible binding for global variable ‘ref_location’ dtMotifMatch: no visible binding for global variable ‘snp_strand’ dtMotifMatch: no visible binding for global variable ‘snp_location’ dtMotifMatch: no visible binding for global variable ‘ref_extra_pwm_left’ dtMotifMatch: no visible binding for global variable ‘ref_extra_pwm_right’ dtMotifMatch: no visible binding for global variable ‘snp_extra_pwm_left’ dtMotifMatch: no visible binding for global variable ‘snp_extra_pwm_right’ dtMotifMatch: no visible binding for global variable ‘snpid’ match_subseq_par: no visible binding for global variable ‘snpid’ match_subseq_par: no visible binding for global variable ‘motif’ match_subseq_par: no visible binding for global variable ‘snpbase’ match_subseq_par: no visible binding for global variable ‘ref_strand’ match_subseq_par: no visible binding for global variable ‘ref_match_seq’ match_subseq_par: no visible binding for global variable ‘ref_seq’ match_subseq_par: no visible binding for global variable ‘ref_start’ match_subseq_par: no visible binding for global variable ‘ref_end’ match_subseq_par: no visible binding for global variable ‘ref_seq_rev’ match_subseq_par: no visible binding for global variable ‘len_seq’ match_subseq_par: no visible binding for global variable ‘snp_strand’ match_subseq_par: no visible binding for global variable ‘snp_match_seq’ match_subseq_par: no visible binding for global variable ‘snp_seq’ match_subseq_par: no visible binding for global variable ‘snp_start’ match_subseq_par: no visible binding for global variable ‘snp_end’ match_subseq_par: no visible binding for global variable ‘snp_seq_rev’ match_subseq_par: no visible binding for global variable ‘snp_seq_ref_match’ match_subseq_par: no visible binding for global variable ‘ref_seq_snp_match’ match_subseq_par: no visible binding for global variable ‘motif_len’ match_subseq_par: no visible binding for global variable ‘log_lik_ref’ match_subseq_par: no visible binding for global variable ‘log_lik_snp’ match_subseq_par: no visible binding for global variable ‘log_lik_ratio’ match_subseq_par: no visible binding for global variable ‘log_enhance_odds’ match_subseq_par: no visible binding for global variable ‘log_reduce_odds’ match_subseq_par: no visible binding for global variable ‘IUPAC’ motif_score_par: no visible binding for global variable ‘motif’ motif_score_par: no visible binding for global variable ‘snpbase’ plotMotifMatch: no visible global function definition for ‘is’ results_motif_par: no visible binding for global variable ‘p.value’ Undefined global functions or variables: IUPAC IUPAC_CODE_MAP is len_seq log_enhance_odds log_lik_ratio log_lik_ref log_lik_snp log_reduce_odds motif motif_len p.value pval_cond_ref pval_cond_snp pval_diff pval_rank pval_ref pval_snp ref_aug_match_seq_forward ref_aug_match_seq_reverse ref_end ref_extra_pwm_left ref_extra_pwm_right ref_location ref_match_seq ref_seq ref_seq_rev ref_seq_snp_match ref_start ref_strand snp_aug_match_seq_forward snp_aug_match_seq_reverse snp_end snp_extra_pwm_left snp_extra_pwm_right snp_location snp_match_seq snp_ref_end snp_ref_length snp_ref_start snp_seq snp_seq_ref_match snp_seq_rev snp_start snp_strand snpbase snpid Consider adding importFrom("methods", "is") to your NAMESPACE file (and ensure that your DESCRIPTION Imports field contains 'methods'). * checking Rd files ... OK * checking Rd metadata ... OK * checking Rd cross-references ... OK * checking for missing documentation entries ... OK * checking for code/documentation mismatches ... OK * checking Rd \usage sections ... OK * checking Rd contents ... OK * checking for unstated dependencies in examples ... OK * checking contents of ‘data’ directory ... OK * checking data for non-ASCII characters ... OK * checking data for ASCII and uncompressed saves ... OK * checking line endings in C/C++/Fortran sources/headers ... OK * checking compiled code ... NOTE Note: information on .o files is not available File ‘/Library/Frameworks/R.framework/Versions/4.3-arm64/Resources/library/atSNP/libs/atSNP.so’: Found ‘_printf’, possibly from ‘printf’ (C) Compiled code should not call entry points which might terminate R nor write to stdout/stderr instead of to the console, nor use Fortran I/O nor system RNGs nor [v]sprintf. The detected symbols are linked into the code but might come from libraries and not actually be called. See ‘Writing portable packages’ in the ‘Writing R Extensions’ manual. * checking files in ‘vignettes’ ... OK * checking examples ... OK * checking for unstated dependencies in ‘tests’ ... OK * checking tests ... Running ‘test.R’ Running ‘test_change.R’ Running ‘test_diff.R’ Running ‘test_is.R’ OK * checking for unstated dependencies in vignettes ... OK * checking package vignettes in ‘inst/doc’ ... OK * checking running R code from vignettes ... SKIPPED * checking re-building of vignette outputs ... SKIPPED * checking PDF version of manual ... OK * DONE Status: 3 NOTEs See ‘/Users/biocbuild/bbs-3.17-bioc-mac-arm64/meat/atSNP.Rcheck/00check.log’ for details.
atSNP.Rcheck/00install.out
############################################################################## ############################################################################## ### ### Running command: ### ### /Library/Frameworks/R.framework/Resources/bin/R CMD INSTALL atSNP ### ############################################################################## ############################################################################## * installing to library ‘/Library/Frameworks/R.framework/Versions/4.3-arm64/Resources/library’ * installing *source* package ‘atSNP’ ... ** using staged installation ** libs using C++ compiler: ‘Apple clang version 14.0.0 (clang-1400.0.29.202)’ using SDK: ‘MacOSX11.3.sdk’ clang++ -arch arm64 -std=gnu++17 -I"/Library/Frameworks/R.framework/Resources/include" -DNDEBUG -I'/Library/Frameworks/R.framework/Versions/4.3-arm64/Resources/library/Rcpp/include' -I/opt/R/arm64/include -fPIC -falign-functions=64 -Wall -g -O2 -c ImportanceSample.cpp -o ImportanceSample.o ImportanceSample.cpp:78:9: warning: variable 'mean_adj_score' set but not used [-Wunused-but-set-variable] double mean_adj_score = 0; ^ ImportanceSample.cpp:77:9: warning: variable 'mean_sample' set but not used [-Wunused-but-set-variable] double mean_sample = 0; ^ 2 warnings generated. clang++ -arch arm64 -std=gnu++17 -I"/Library/Frameworks/R.framework/Resources/include" -DNDEBUG -I'/Library/Frameworks/R.framework/Versions/4.3-arm64/Resources/library/Rcpp/include' -I/opt/R/arm64/include -fPIC -falign-functions=64 -Wall -g -O2 -c ImportanceSampleChange.cpp -o ImportanceSampleChange.o ImportanceSampleChange.cpp:49:9: warning: variable 'mean_score' set but not used [-Wunused-but-set-variable] double mean_score = 0; ^ 1 warning generated. clang++ -arch arm64 -std=gnu++17 -I"/Library/Frameworks/R.framework/Resources/include" -DNDEBUG -I'/Library/Frameworks/R.framework/Versions/4.3-arm64/Resources/library/Rcpp/include' -I/opt/R/arm64/include -fPIC -falign-functions=64 -Wall -g -O2 -c ImportanceSampleDiff.cpp -o ImportanceSampleDiff.o ImportanceSampleDiff.cpp:85:9: warning: variable 'mean_diff' set but not used [-Wunused-but-set-variable] double mean_diff = 0; ^ 1 warning generated. clang++ -arch arm64 -std=gnu++17 -I"/Library/Frameworks/R.framework/Resources/include" -DNDEBUG -I'/Library/Frameworks/R.framework/Versions/4.3-arm64/Resources/library/Rcpp/include' -I/opt/R/arm64/include -fPIC -falign-functions=64 -Wall -g -O2 -c MotifScore.cpp -o MotifScore.o clang++ -arch arm64 -std=gnu++17 -dynamiclib -Wl,-headerpad_max_install_names -undefined dynamic_lookup -single_module -multiply_defined suppress -L/Library/Frameworks/R.framework/Resources/lib -L/opt/R/arm64/lib -o atSNP.so ImportanceSample.o ImportanceSampleChange.o ImportanceSampleDiff.o MotifScore.o -F/Library/Frameworks/R.framework/.. -framework R -Wl,-framework -Wl,CoreFoundation installing to /Library/Frameworks/R.framework/Versions/4.3-arm64/Resources/library/00LOCK-atSNP/00new/atSNP/libs ** R ** data ** inst ** byte-compile and prepare package for lazy loading ** help *** installing help indices ** building package indices ** installing vignettes ** testing if installed package can be loaded from temporary location ** checking absolute paths in shared objects and dynamic libraries ** testing if installed package can be loaded from final location ** testing if installed package keeps a record of temporary installation path * DONE (atSNP)
atSNP.Rcheck/tests/test.Rout
R version 4.3.1 (2023-06-16) -- "Beagle Scouts" Copyright (C) 2023 The R Foundation for Statistical Computing Platform: aarch64-apple-darwin20 (64-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > library(atSNP) > library(BiocParallel) > library(testthat) > > ## process the data > data(example) > > motif_scores <- ComputeMotifScore(motif_library, snpInfo, ncores = 1) > > motif_scores <- MatchSubsequence(motif_scores$snp.tbl, motif_scores$motif.scores, ncores = 1, motif.lib = motif_library) > > motif_scores[which(motif_scores$snpid == "rs7412" & motif_scores$motif == "SIX5_disc1"), ] snpid motif 4 rs7412 SIX5_disc1 ref_seq 4 CTCCTCCGCGATGCCGATGACCTGCAGAAGCGCCTGGCAGTGTACCAGGCCGGGGCCCGCG snp_seq motif_len 4 CTCCTCCGCGATGCCGATGACCTGCAGAAGTGCCTGGCAGTGTACCAGGCCGGGGCCCGCG 10 ref_start ref_end ref_strand snp_start snp_end snp_strand log_lik_ref 4 29 38 - 22 31 + -42.60672 log_lik_snp log_lik_ratio log_enhance_odds log_reduce_odds IUPAC 4 -38.4083 -4.198418 23.013 -2.917768 GARWTGTAGT ref_match_seq snp_match_seq ref_seq_snp_match snp_seq_ref_match snpbase 4 GCCAGGCGCT CTGCAGAAGT CTGCAGAAGC GCCAGGCACT T > > len_seq <- sapply(motif_scores$ref_seq, nchar) > snp_pos <- as.integer(len_seq / 2) + 1 > > i <- which(motif_scores$snpid == "rs7412" & motif_scores$motif == "SIX5_disc1") > > test_that("Error: reference bases are not the same as the sequence matrix.", { + expect_equal(sum(snpInfo$sequence_matrix[31, ] != snpInfo$ref_base), 0) + expect_equal(sum(snpInfo$sequence_matrix[31, ] == snpInfo$snp_base), 0) + }) Test passed 🥇 > > test_that("Error: log_lik_ratio is not correct.", { + expect_equal(motif_scores$log_lik_ref - motif_scores$log_lik_snp, motif_scores$log_lik_ratio) + }) Test passed 🥳 > > test_that("Error: log likelihoods are not correct.", { + + log_lik <- sapply(seq(nrow(motif_scores)), + function(i) { + motif_mat <- motif_library[[motif_scores$motif[i]]] + colind<-which(snpInfo$snpids==motif_scores$snpid[i]) + bases <- snpInfo$sequence_matrix[motif_scores$ref_start[i]:motif_scores$ref_end[i], colind] + if(motif_scores$ref_strand[i] == "-") + bases <- 5 - rev(bases) + log(prod( + motif_mat[cbind(seq(nrow(motif_mat)), + bases)])) + }) + + expect_equal(log_lik, motif_scores$log_lik_ref) + + snp_mat <- snpInfo$sequence_matrix + snp_mat[cbind(snp_pos, seq(ncol(snp_mat)))] <- snpInfo$snp_base + log_lik <- sapply(seq(nrow(motif_scores)), + function(i) { + motif_mat <- motif_library[[motif_scores$motif[i]]] + colind<-which(snpInfo$snpids==motif_scores$snpid[i]) + bases <- snp_mat[motif_scores$snp_start[i]:motif_scores$snp_end[i], colind] + if(motif_scores$snp_strand[i] == "-") + bases <- 5 - rev(bases) + log(prod( + motif_mat[cbind(seq(nrow(motif_mat)), + bases)])) + }) + + expect_equal(log_lik, motif_scores$log_lik_snp) + }) Test passed 😀 > > test_that("Error: log_enhance_odds not correct.", { + + len_seq <- sapply(motif_scores$ref_seq, nchar) + snp_pos <- as.integer(len_seq / 2) + 1 + + ## log odds for reduction in binding affinity + + pos_in_pwm <- snp_pos - motif_scores$ref_start + 1 + neg_ids <- which(motif_scores$ref_strand == "-") + pos_in_pwm[neg_ids] <- motif_scores$ref_end[neg_ids]- snp_pos[neg_ids] + 1 + snp_base <- sapply(substr(motif_scores$snp_seq, snp_pos, snp_pos), function(x) which(c("A", "C", "G", "T") == x)) + ref_base <- sapply(substr(motif_scores$ref_seq, snp_pos, snp_pos), function(x) which(c("A", "C", "G", "T") == x)) + snp_base[neg_ids] <- 5 - snp_base[neg_ids] + ref_base[neg_ids] <- 5 - ref_base[neg_ids] + my_log_reduce_odds <- sapply(seq(nrow(motif_scores)), + function(i) + log(motif_library[[motif_scores$motif[i]]][pos_in_pwm[i], ref_base[i]]) - + log(motif_library[[motif_scores$motif[i]]][pos_in_pwm[i], snp_base[i]]) + ) + + expect_equal(my_log_reduce_odds, motif_scores$log_reduce_odds) + + ## log odds in enhancing binding affinity + + pos_in_pwm <- snp_pos - motif_scores$snp_start + 1 + neg_ids <- which(motif_scores$snp_strand == "-") + pos_in_pwm[neg_ids] <- motif_scores$snp_end[neg_ids]- snp_pos[neg_ids] + 1 + snp_base <- sapply(substr(motif_scores$snp_seq, snp_pos, snp_pos), function(x) which(c("A", "C", "G", "T") == x)) + ref_base <- sapply(substr(motif_scores$ref_seq, snp_pos, snp_pos), function(x) which(c("A", "C", "G", "T") == x)) + snp_base[neg_ids] <- 5 - snp_base[neg_ids] + ref_base[neg_ids] <- 5 - ref_base[neg_ids] + my_log_enhance_odds <- sapply(seq(nrow(motif_scores)), + function(i) + log(motif_library[[motif_scores$motif[i]]][pos_in_pwm[i], snp_base[i]]) - + log(motif_library[[motif_scores$motif[i]]][pos_in_pwm[i], ref_base[i]]) + ) + + expect_equal(my_log_enhance_odds, motif_scores$log_enhance_odds) + + + }) Test passed 🎉 > > test_that("Error: the maximum log likelihood computation is not correct.", { + + snp_mat <- snpInfo$sequence_matrix + snp_mat[cbind(snp_pos, seq(ncol(snp_mat)))] <- snpInfo$snp_base + + .findMaxLog <- function(seq_vec, pwm) { + snp_pos <- as.integer(length(seq_vec) / 2) + 1 + start_pos <- snp_pos - nrow(pwm) + 1 + end_pos <- snp_pos + rev_seq <- 5 - rev(seq_vec) + + maxLogProb <- -Inf + for(i in start_pos : end_pos) { + LogProb <- log(prod(pwm[cbind(seq(nrow(pwm)), + seq_vec[i - 1 + seq(nrow(pwm))])])) + if(LogProb > maxLogProb) + maxLogProb <- LogProb + } + for(i in start_pos : end_pos) { + LogProb <- log(prod(pwm[cbind(seq(nrow(pwm)), + rev_seq[i - 1 + seq(nrow(pwm))])])) + if(LogProb > maxLogProb) + maxLogProb <- LogProb + } + return(maxLogProb) + } + + ## find the maximum log likelihood on the reference sequence + my_log_lik_ref <- sapply(seq(nrow(motif_scores)), + function(x) { + colind<-which(snpInfo$snpids==motif_scores$snpid[x]) + seq_vec<- snpInfo$sequence_matrix[, colind] + pwm <- motif_library[[motif_scores$motif[x]]] + return(.findMaxLog(seq_vec, pwm)) + }) + + ## find the maximum log likelihood on the SNP sequence + + my_log_lik_snp <- sapply(seq(nrow(motif_scores)), + function(x) { + colind<-which(snpInfo$snpids==motif_scores$snpid[x]) #ADDED + seq_vec<- snp_mat[, colind] + pwm <- motif_library[[motif_scores$motif[x]]] + return(.findMaxLog(seq_vec, pwm)) + }) + + expect_equal(my_log_lik_ref, motif_scores$log_lik_ref) + expect_equal(my_log_lik_snp, motif_scores$log_lik_snp) + + }) Test passed 😸 > > proc.time() user system elapsed 16.119 0.820 26.104
atSNP.Rcheck/tests/test_change.Rout
R version 4.3.1 (2023-06-16) -- "Beagle Scouts" Copyright (C) 2023 The R Foundation for Statistical Computing Platform: aarch64-apple-darwin20 (64-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > library(atSNP) > library(BiocParallel) > library(testthat) > data(example) > > trans_mat <- matrix(rep(snpInfo$prior, each = 4), nrow = 4) > test_pwm <- motif_library$SIX5_disc1 > scores <- as.matrix(motif_scores$motif.scores[3:4, 4:5]) > score_diff <- abs(scores[,2]-scores[,1]) > > pval_a <- .Call("test_p_value", test_pwm, snpInfo$prior, snpInfo$transition, scores, 0.15, 100) > pval_ratio <- abs(log(pval_a[seq(nrow(scores)),1]) - log(pval_a[seq(nrow(scores)) + nrow(scores), 1])) > > test_score <- test_pwm > for(i in seq(nrow(test_score))) { + for(j in seq(ncol(test_score))) { + test_score[i, j] <- exp(mean(log(test_pwm[i, j] / test_pwm[i, -j]))) + } + } > > adj_mat <- test_pwm + 0.25 > motif_len <- nrow(test_pwm) > > ## these are functions for this test only > drawonesample <- function(theta) { + prob_start <- rev(rowSums(test_score ^ theta) / rowSums(adj_mat)) + id <- sample(seq(motif_len), 1, prob = prob_start) + sample <- sample(1:4, 2 * motif_len - 1, replace = TRUE, prob = snpInfo$prior) + delta <- adj_mat + delta[motif_len - id + 1, ] <- test_score[motif_len - id + 1, ] ^ theta + sample[id - 1 + seq(motif_len)] <- apply(delta, 1, function(x) sample(seq(4), 1, prob = x)) + ## compute weight + sc <- 0 + for(s in seq(motif_len)) { + delta <- adj_mat + delta[motif_len + 1 - s, ] <- test_score[motif_len + 1 - s, ] ^ theta + sc <- sc + prod(delta[cbind(seq(motif_len), sample[s - 1 + seq(motif_len)])]) / + prod(snpInfo$prior[sample[s - 1 + seq(motif_len)]]) + } + sample <- c(sample, id, sc) + return(sample) + } > jointprob <- function(x) prod(test_pwm[cbind(seq(motif_len), x)]) > maxjointprob <- function(x) { + maxp <- -Inf + p <- -Inf + for(i in 1:motif_len) { + p <- jointprob(x[i:(i+motif_len - 1)]) + if(p > maxp) + maxp <- p + } + for(i in 1:motif_len) { + p <- jointprob(5 - x[(i+motif_len - 1):i]) + if(p > maxp) + maxp <- p + } + return(maxp) + } > get_freq <- function(sample) { + emp_freq <- matrix(0, nrow = 2 * motif_len - 1, ncol = 4) + for(i in seq(2 * motif_len - 1)) { + for(j in seq(4)) { + emp_freq[i, j] <- sum(sample[i, ] == j - 1) + } + } + emp_freq <- emp_freq / rowSums(emp_freq) + return(emp_freq) + } > > test_that("Error: quantile function computing are not equivalent.", { + for(p in c(0.01, 0.1, 0.5, 0.9, 0.99) ) { + delta <- .Call("test_find_percentile_change", score_diff, p, package = "atSNP") + delta.r <- as.double(sort(abs(scores[,2]-scores[,1]))[ceiling((1 - p) * (nrow(scores)))]) + expect_equal(delta, delta.r) + } + }) Test passed 🎉 > > test_that("Error: the scores for samples are not equivalent.", { + p <- 0.1 + delta <- .Call("test_find_percentile_change", score_diff, p, package = "atSNP") + theta <- .Call("test_find_theta_change", test_score, adj_mat, delta, package = "atSNP") + ## Use R code to generate a random sample + for(i in seq(10)) { + sample <- drawonesample(theta) + sample_score <- .Call("test_compute_sample_score_change", test_pwm, test_score, adj_mat, sample[seq(2 * motif_len - 1)] - 1, snpInfo$prior, trans_mat, sample[2 * motif_len] - 1, theta, package = "atSNP") + expect_equal(sample[2 * motif_len + 1], sample_score[1]) + sample1 <- sample2 <- sample3 <- sample + sample1[motif_len] <- seq(4)[-sample[motif_len]][1] + sample2[motif_len] <- seq(4)[-sample[motif_len]][2] + sample3[motif_len] <- seq(4)[-sample[motif_len]][3] + sample_score_r <- log(maxjointprob(sample[seq(2 * motif_len - 1)])) - + log(c(maxjointprob(sample1[seq(2 * motif_len - 1)]), + maxjointprob(sample2[seq(2 * motif_len - 1)]), + maxjointprob(sample3[seq(2 * motif_len - 1)]))) + expect_equal(sample_score_r, sample_score[2:4]) + } + + ## Use C code to generate a random sample + for(i in seq(10)) { + sample <- .Call("test_importance_sample_change", test_score, snpInfo$prior, trans_mat, test_pwm, theta, package = "atSNP") + start_pos <- sample[2 * motif_len] + 1 + adj_score <- 0 + for(s in seq_len(motif_len)) { + adj_s <- sum(log(adj_mat[cbind(seq(motif_len), sample[s - 1 + seq(motif_len)] + 1)]) - + log(snpInfo$prior[sample[s - 1 + seq(motif_len)] + 1])) + adj_s <- adj_s + theta * log(test_score[motif_len + 1 - s, sample[motif_len] + 1]) - + log(adj_mat[motif_len + 1 - s, sample[motif_len] + 1]) + adj_score <- adj_score + exp(adj_s) + } + sample_score <- .Call("test_compute_sample_score_change", test_pwm, test_score, adj_mat, sample[seq(2 * motif_len - 1)], snpInfo$prior, trans_mat, sample[2 * motif_len], theta, package = "atSNP") + expect_equal(adj_score, sample_score[1]) + } + }) Test passed 🥳 > > test_that("Error: compute the normalizing constant.", { + ## parameters + for(p in seq(9) / 10) { + delta <- .Call("test_find_percentile_change", score_diff, p, package = "atSNP") + theta <- .Call("test_find_theta_change", test_score, adj_mat, delta, package = "atSNP") + const <- .Call("test_func_delta_change", test_score, adj_mat, theta, package = "atSNP") + ## in R + adj_sum <- rowSums(adj_mat) + wei_sum <- rowSums(test_score ^ theta) + const.r <- prod(adj_sum) * sum(wei_sum / adj_sum) + expect_equal(const, const.r) + } + }) Test passed 😸 > > test_that("Error: sample distributions are not expected.", { + ## parameters + p <- 0.1 + delta <- .Call("test_find_percentile_change", score_diff, p, package = "atSNP") + theta <- .Call("test_find_theta_change", test_score, adj_mat, delta, package = "atSNP") + prob_start <- rev(rowSums(test_score ^ theta) / rowSums(adj_mat)) + ## construct the delta matrix + delta <- matrix(1, nrow = 4 * motif_len, ncol = 2 * motif_len - 1) + for(pos in seq(motif_len)) { + delta[seq(4) + 4 * (pos - 1), ] <- snpInfo$prior + delta[seq(4) + 4 * (pos - 1), pos - 1 + seq(motif_len)] <- t(test_pwm) + delta[seq(4) + 4 * (pos - 1), motif_len] <- test_score[motif_len + 1 - pos, ] ^ theta + delta[seq(4) + 4 * (pos - 1), ] <- delta[seq(4) + 4 * (pos - 1),] / rep(colSums(delta[seq(4) + 4 * (pos - 1), ]), each = 4) + } + target_freq <- matrix(0, nrow = 4, ncol = 2 * motif_len - 1) + for(pos in seq(motif_len)) { + target_freq <- target_freq + delta[seq(4) + 4 * (pos - 1), ] * prob_start[pos] + } + target_freq <- t(target_freq) + target_freq <- target_freq / rowSums(target_freq) + + results_i <- function(i) { + ## generate 100 samples + sample1 <- sapply(seq(100), function(x) + .Call("test_importance_sample_change", + adj_mat, snpInfo$prior, trans_mat, test_score, theta, package = "atSNP")) + emp_freq1 <- get_freq(sample1) + sample2 <- sapply(rep(theta, 100), drawonesample) + emp_freq2 <- get_freq(sample2 - 1) + ## print(rbind(emp_freq1[10, ], emp_freq2[10, ], target_freq[10, ])) + max(abs(emp_freq1 - target_freq)) > max(abs(emp_freq2 - target_freq)) + } + + if(Sys.info()[["sysname"]] == "Windows"){ + snow <- SnowParam(workers = 1, type = "SOCK") + results<-bpmapply(results_i, seq(20), BPPARAM = snow,SIMPLIFY = FALSE) + }else{ + results<-bpmapply(results_i, seq(20), BPPARAM = MulticoreParam(workers = 1), + SIMPLIFY = FALSE) + } + + print(sum(unlist(results))) + print(pbinom(sum(unlist(results)), size = 20, prob = 0.5)) + }) [1] 9 [1] 0.4119015 ── Skip (???): Error: sample distributions are not expected. ─────────────────── Reason: empty test > > test_that("Error: the chosen pvalues should have the smaller variance.", { + .structure_diff <- function(pval_mat) { + id <- apply(pval_mat[, c(2, 4)], 1, which.min) + return(cbind(pval_mat[, c(1, 3)][cbind(seq_along(id), id)], + pval_mat[, c(2, 4)][cbind(seq_along(id), id)])) + } + for(p in c(0.05, 0.1, 0.2, 0.5)) { + p_values <- .Call("test_p_value_change", test_pwm, test_score, adj_mat, snpInfo$prior, snpInfo$transition, score_diff, pval_ratio, quantile(score_diff, 1 - p), 100, package = "atSNP")$score + p_values_s <- .structure_diff(p_values) + expect_equal(p_values_s[, 2], apply(p_values[, c(2, 4)], 1, min)) + } + }) Test passed 😀 > > proc.time() user system elapsed 18.690 0.889 28.756
atSNP.Rcheck/tests/test_diff.Rout
R version 4.3.1 (2023-06-16) -- "Beagle Scouts" Copyright (C) 2023 The R Foundation for Statistical Computing Platform: aarch64-apple-darwin20 (64-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > library(atSNP) > library(BiocParallel) > library(testthat) > data(example) > > trans_mat <- matrix(rep(snpInfo$prior, each = 4), nrow = 4) > test_pwm <- motif_library$SIX5_disc1 > scores <- as.matrix(motif_scores$motif.scores[3:4, 4:5]) > score_diff <- abs(scores[,2]-scores[,1]) > > test_score <- test_pwm > for(i in seq(nrow(test_score))) { + for(j in seq(ncol(test_score))) { + test_score[i, j] <- exp(mean(log(test_pwm[i, j] / test_pwm[i, -j]))) + } + } > > adj_mat <- test_pwm + rowMeans(test_pwm) > motif_len <- nrow(test_pwm) > > ## these are functions for this test only > drawonesample <- function(theta) { + prob_start <- sapply(seq(motif_len), + function(j) + sum(snpInfo$prior * test_score[motif_len + 1 - j, ] ^ theta * + adj_mat[motif_len + 1 - j, ]) / + sum(snpInfo$prior * adj_mat[motif_len + 1 - j, ]) + ) + id <- sample(seq(motif_len), 1, prob = prob_start) + sample <- sample(1:4, 2 * motif_len - 1, replace = TRUE, prob = snpInfo$prior) + delta <- adj_mat + delta[motif_len + 1 - id, ] <- delta[motif_len + 1 - id, ] * test_score[motif_len + 1 - id, ] ^ theta + sample[id - 1 + seq(motif_len)] <- apply(delta, 1, function(x) + sample(seq(4), 1, prob = x * snpInfo$prior)) + sc <- 0 + for(s in seq(motif_len)) { + delta <- adj_mat + delta[motif_len + 1 - s, ] <- delta[motif_len + 1 - s, ] * test_score[motif_len + 1 - s, ] ^ theta + sc <- sc + prod(delta[cbind(seq(motif_len), sample[s - 1 + seq(motif_len)])]) + } + sample <- c(sample, id, sc) + return(sample) + } > jointprob <- function(x) prod(test_pwm[cbind(seq(motif_len), x)]) > maxjointprob <- function(x) { + maxp <- -Inf + p <- -Inf + for(i in 1:motif_len) { + p <- jointprob(x[i:(i+motif_len - 1)]) + if(p > maxp) + maxp <- p + } + for(i in 1:motif_len) { + p <- jointprob(5 - x[(i+motif_len - 1):i]) + if(p > maxp) + maxp <- p + } + return(maxp) + } > get_freq <- function(sample) { + emp_freq <- matrix(0, nrow = 2 * motif_len - 1, ncol = 4) + for(i in seq(2 * motif_len - 1)) { + for(j in seq(4)) { + emp_freq[i, j] <- sum(sample[i, ] == j - 1) + } + } + emp_freq <- emp_freq / rowSums(emp_freq) + return(emp_freq) + } > > test_that("Error: quantile function computing are not equivalent.", { + for(p in c(0.01, 0.1, 0.5, 0.9, 0.99)) { + delta <- .Call("test_find_percentile_diff", score_diff, p, package = "atSNP") + delta.r <- as.double(sort(abs(scores[,2]-scores[,1]))[ceiling((1 - p) * (nrow(scores)))]) + expect_equal(delta, delta.r) + } + }) Test passed 🥇 > > test_that("Error: the scores for samples are not equivalent.", { + p <- 0.1 + delta <- .Call("test_find_percentile_diff", score_diff, p, package = "atSNP") + theta <- .Call("test_find_theta_diff", test_score, adj_mat, snpInfo$prior, snpInfo$transition, delta, package = "atSNP") + ## Use R code to generate a random sample + for(i in seq(10)) { + sample <- drawonesample(theta) + sample_score <- .Call("test_compute_sample_score_diff", test_pwm, test_score, adj_mat, sample[seq(2 * motif_len - 1)] - 1, sample[2 * motif_len] - 1, theta, package = "atSNP") + expect_equal(sample[2 * motif_len + 1], sample_score[1]) + sample1 <- sample2 <- sample3 <- sample + sample1[motif_len] <- seq(4)[-sample[motif_len]][1] + sample2[motif_len] <- seq(4)[-sample[motif_len]][2] + sample3[motif_len] <- seq(4)[-sample[motif_len]][3] + sample_score_r <- log(maxjointprob(sample[seq(2 * motif_len - 1)])) - + log(c(maxjointprob(sample1[seq(2 * motif_len - 1)]), + maxjointprob(sample2[seq(2 * motif_len - 1)]), + maxjointprob(sample3[seq(2 * motif_len - 1)]))) + expect_equal(sample_score_r, sample_score[-1]) + } + + ## Use C code to generate a random sample + delta <- matrix(1, nrow = 4 * motif_len, ncol = 2 * motif_len - 1) + for(pos in seq(motif_len)) { + for(j in (pos + motif_len - 1) : 1) { + if(j < pos + motif_len - 1) { + delta[4 * (pos - 1) + seq(4), j] <- sum(snpInfo$prior * delta[4 * (pos - 1) + seq(4), j + 1]) + } + if(j >= pos) { + delta[4 * (pos - 1) + seq(4), j] <- delta[4 * (pos - 1) + seq(4), j] * adj_mat[j - pos + 1, ] + } + if(j == motif_len) { + delta[4 * (pos - 1) + seq(4), j] <- delta[4 * (pos - 1) + seq(4), j] * test_score[j - pos + 1, ] ^ theta + } + } + } + for(i in seq(10)) { + sample <- .Call("test_importance_sample_diff", delta, snpInfo$prior, trans_mat, test_pwm, theta, package = "atSNP") + start_pos <- sample[2 * motif_len] + 1 + adj_score <- 0 + for(s in seq_len(motif_len)) { + adj_s <- sum(log(adj_mat[cbind(seq(motif_len), sample[s - 1 + seq(motif_len)] + 1)])) + adj_s <- adj_s + theta * log(test_score[motif_len + 1 - s, sample[motif_len] + 1]) + adj_score <- adj_score + exp(adj_s) + } + sample_score <- .Call("test_compute_sample_score_diff", test_pwm, test_score, adj_mat, sample[seq(2 * motif_len - 1)], sample[2 * motif_len], theta, package = "atSNP") + expect_equal(adj_score, sample_score[1]) + } + }) Test passed 😀 > > test_that("Error: compute the normalizing constant.", { + + ## parameters + p <- 0.1 + delta <- .Call("test_find_percentile_diff", score_diff, p, package = "atSNP") + theta <- .Call("test_find_theta_diff", test_score, adj_mat, snpInfo$prior, snpInfo$transition, delta, package = "atSNP") + + ## + const <- .Call("test_func_delta_diff", test_score, adj_mat, snpInfo$prior, trans_mat, theta, package = "atSNP") + + prob_start <- sapply(seq(motif_len), + function(j) + sum(snpInfo$prior * test_score[motif_len + 1 - j, ] ^ theta * + adj_mat[motif_len + 1 - j, ]) / + sum(snpInfo$prior * adj_mat[motif_len + 1 - j, ]) + ) + + const.r <- prod(colSums(snpInfo$prior * t(adj_mat))) * sum(prob_start) + expect_equal(const, const.r) + }) Test passed 😸 > > test_that("Error: sample distributions are not expected.", { + + ## parameters + p <- 0.1 + delta <- .Call("test_find_percentile_diff", score_diff, p, package = "atSNP") + theta <- .Call("test_find_theta_diff", test_score, adj_mat, snpInfo$prior, snpInfo$transition, delta, package = "atSNP") + + ## construct the delta matrix + delta <- matrix(1, nrow = 4 * motif_len, ncol = 2 * motif_len - 1) + for(pos in seq(motif_len)) { + for(j in (pos + motif_len - 1) : 1) { + if(j < pos + motif_len - 1) { + delta[4 * (pos - 1) + seq(4), j] <- sum(snpInfo$prior * delta[4 * (pos - 1) + seq(4), j + 1]) + } + if(j >= pos) { + delta[4 * (pos - 1) + seq(4), j] <- delta[4 * (pos - 1) + seq(4), j] * adj_mat[j - pos + 1, ] + } + if(j == motif_len) { + delta[4 * (pos - 1) + seq(4), j] <- delta[4 * (pos - 1) + seq(4), j] * test_score[j - pos + 1, ] ^ theta + } + } + } + + target_freq <- matrix(0, nrow = 4, ncol = 2 * motif_len - 1) + + mat <- snpInfo$prior * matrix(delta[, 1], nrow = 4) + wei <- colSums(mat) + for(j in seq(2 * motif_len - 1)) { + for(pos in seq(motif_len)) { + tmp <- delta[seq(4) + 4 * (pos - 1), j] * snpInfo$prior + target_freq[, j] <- target_freq[, j] + tmp / sum(tmp) * wei[pos] + } + } + target_freq <- t(target_freq) + target_freq <- target_freq / rowSums(target_freq) + + results_i <- function(i) { + ## generate 100 samples + sample1 <- sapply(seq(100), function(x) + .Call("test_importance_sample_diff", + delta, snpInfo$prior, trans_mat, test_score, theta, package = "atSNP")) + emp_freq1 <- get_freq(sample1) + + sample2 <- sapply(rep(theta, 100), drawonesample) + emp_freq2 <- get_freq(sample2 - 1) + + ## print(rbind(emp_freq1[10, ], emp_freq2[10, ], target_freq[10, ])) + max(abs(emp_freq1 - target_freq)) > max(abs(emp_freq2 - target_freq)) + } + + if(Sys.info()[["sysname"]] == "Windows"){ + snow <- SnowParam(workers = 1, type = "SOCK") + results<-bpmapply(results_i, seq(20), BPPARAM = snow,SIMPLIFY = FALSE) + }else{ + results<-bpmapply(results_i, seq(20), BPPARAM = MulticoreParam(workers = 1), + SIMPLIFY = FALSE) + } + + print(sum(unlist(results))) + + print(pbinom(sum(unlist(results)), size = 20, prob = 0.5)) + + }) [1] 11 [1] 0.7482777 ── Skip (???): Error: sample distributions are not expected. ─────────────────── Reason: empty test > > test_that("Error: the chosen pvalues should have the smaller variance.", { + + .structure_diff <- function(pval_mat) { + id <- apply(pval_mat[, c(2, 4)], 1, which.min) + return(cbind(pval_mat[, c(1, 3)][cbind(seq_along(id), id)], + pval_mat[, c(2, 4)][cbind(seq_along(id), id)])) + } + + for(p in c(0.05, 0.1, 0.2, 0.5)) { + p_values <- .Call("test_p_value_diff", test_pwm, test_score, adj_mat, snpInfo$prior, snpInfo$transition, score_diff, quantile(score_diff, 1 - p), 100, package = "atSNP") + p_values_s <- .structure_diff(p_values) + expect_equal(p_values_s[, 2], apply(p_values[, c(2, 4)], 1, min)) + } + }) Test passed 🎉 > > proc.time() user system elapsed 18.704 0.852 29.106
atSNP.Rcheck/tests/test_is.Rout
R version 4.3.1 (2023-06-16) -- "Beagle Scouts" Copyright (C) 2023 The R Foundation for Statistical Computing Platform: aarch64-apple-darwin20 (64-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > library(atSNP) > library(BiocParallel) > library(testthat) > data(example) > > trans_mat <- matrix(rep(snpInfo$prior, each = 4), nrow = 4) > test_pwm <- motif_library$SIX5_disc1 > scores <- as.matrix(motif_scores$motif.scores[3:4, 4:5]) > > motif_len <- nrow(test_pwm) > > ## these are functions for this test only > drawonesample <- function(theta) { + delta <- snpInfo$prior * t(test_pwm ^ theta) + delta <- delta / rep(colSums(delta), each = 4) + sample <- sample(1:4, 2 * motif_len - 1, replace = TRUE, prob = snpInfo$prior) + id <- sample(seq(motif_len), 1) + sample[id : (id + motif_len - 1)] <- apply(delta, 2, function(x) sample(1:4, 1, prob = x)) + sc <- s_cond <- 0 + for(s in seq(motif_len)) { + sc <- sc + prod(test_pwm[cbind(seq(motif_len), + sample[s : (s + motif_len - 1)])]) ^ theta + } + s_cond <- prod(test_pwm[cbind(seq(motif_len), + sample[id : (id + motif_len - 1)])]) ^ theta + sample <- c(sample, id, sc, s_cond) + return(sample) + } > jointprob <- function(x) prod(test_pwm[cbind(seq(motif_len), x)]) > maxjointprob <- function(x) { + maxp <- -Inf + p <- -Inf + for(i in 1:motif_len) { + p <- jointprob(x[i:(i+motif_len - 1)]) + if(p > maxp) + maxp <- p + } + for(i in 1:motif_len) { + p <- jointprob(5 - x[(i+motif_len - 1):i]) + if(p > maxp) + maxp <- p + } + return(maxp) + } > get_freq <- function(sample) { + ids <- cbind( + rep(sample[motif_len * 2, ], each = motif_len) + seq(motif_len), + rep(seq(100), each = motif_len)) + sample_motif <- matrix(sample[ids], nrow = motif_len) + 1 + emp_freq <- matrix(0, nrow = motif_len, ncol = 4) + for(i in seq(motif_len)) { + for(j in seq(4)) { + emp_freq[i, j] <- sum(sample_motif[i, ] == j) + } + } + emp_freq <- emp_freq / rowSums(emp_freq) + return(emp_freq) + } > > test_that("Error: quantile function computing are not equivalent.", { + for(p in c(0.01, 0.1, 0.5, 0.9, 0.99)) { + delta <- .Call("test_find_percentile", c(scores), p, package = "atSNP") + delta.r <- -sort(-c(scores))[as.integer(p * length(scores)) + 1] + expect_equal(delta, delta.r) + } + }) Test passed 🌈 > > test_that("Error: the scores for samples are not equivalent.", { + p <- 0.01 + delta <- .Call("test_find_percentile", scores, p, package = "atSNP") + theta <- .Call("test_find_theta", test_pwm, snpInfo$prior, snpInfo$transition, delta, package = "atSNP") + ## Use R code to generate a random sample + for(i in seq(10)) { + sample <- drawonesample(theta) + sample_score <- .Call("test_compute_sample_score", test_pwm, sample[seq(2 * motif_len - 1)] - 1, sample[motif_len * 2] - 1, theta, package = "atSNP") + expect_equal(sample[2 * motif_len + 1], sample_score[2]) + expect_equal(sample[2 * motif_len + 2], sample_score[3]) + } + ## Use C code to generate a random sample + for(i in seq(10)) { + delta <- t(test_pwm ^ theta) + delta <- cbind(matrix( + sum(snpInfo$prior * delta[, 1]), + nrow = 4, ncol = motif_len - 1), delta) + sample <- .Call("test_importance_sample", delta, snpInfo$prior, trans_mat, test_pwm, theta, package = "atSNP") + start_pos <- sample[motif_len * 2] + adj_score <- 0 + for(s in seq(motif_len) - 1) { + adj_score <- adj_score + prod(test_pwm[cbind(seq(motif_len), + sample[s + seq(motif_len)] + 1)]) ^ theta + } + adj_score_cond <- prod(test_pwm[cbind(seq(motif_len), sample[start_pos + seq(motif_len)] + 1)]) ^ theta + sample_score <- .Call("test_compute_sample_score", test_pwm, sample[seq(2 * motif_len - 1)], sample[motif_len * 2], theta, package = "atSNP") + expect_equal(adj_score, sample_score[2]) + expect_equal(adj_score_cond, sample_score[3]) + } + }) Test passed 😸 > > test_that("Error: compute the normalizing constant.", { + ## parameters + p <- 0.01 + delta <- .Call("test_find_percentile", scores, p, package = "atSNP") + theta <- .Call("test_find_theta", test_pwm, snpInfo$prior, snpInfo$transition, delta, package = "atSNP") + ## + const <- .Call("test_func_delta", test_pwm, snpInfo$prior, trans_mat, theta, package = "atSNP") + const.r <- prod(colSums(snpInfo$prior * t(test_pwm) ^ theta)) * motif_len + expect_equal(abs(const - const.r) / const < 1e-5, TRUE) + }) Test passed 🥳 > > test_that("Error: sample distributions are not expected.", { + ## parameters + p <- 0.1 + delta <- .Call("test_find_percentile", scores, p, package = "atSNP") + theta <- .Call("test_find_theta", test_pwm, snpInfo$prior, trans_mat, delta, package = "atSNP") + delta <- t(test_pwm ^ theta) + delta <- cbind(matrix( + sum(snpInfo$prior * delta[, 1]), + nrow = 4, ncol = motif_len - 1), delta) + + results_i <- function(i) { + ## generate 100 samples + sample <- sapply(seq(100), function(x) + .Call("test_importance_sample", + delta, snpInfo$prior, trans_mat, test_pwm, theta, package = "atSNP")) + emp_freq1 <- get_freq(sample) + target_freq <- test_pwm ^ theta * snpInfo$prior + target_freq <- target_freq / rowSums(target_freq) + ## generate samples in R + sample <- sapply(rep(theta, 100), drawonesample) + emp_freq2 <- get_freq(sample[seq(2 * motif_len), ] - 1) + max(abs(emp_freq1 - target_freq)) > max(abs(emp_freq2 - target_freq)) + } + + if(Sys.info()[["sysname"]] == "Windows"){ + snow <- SnowParam(workers = 1, type = "SOCK") + results<-bpmapply(results_i, seq(20), BPPARAM = snow,SIMPLIFY = FALSE) + }else{ + results<-bpmapply(results_i, seq(20), BPPARAM = MulticoreParam(workers = 1), + SIMPLIFY = FALSE) + } + + print(sum(unlist(results))) + print(pbinom(sum(unlist(results)), size = 20, prob = 0.5)) + }) [1] 10 [1] 0.5880985 ── Skip (???): Error: sample distributions are not expected. ─────────────────── Reason: empty test > > test_that("Error: the chosen pvalues should have the smaller variance.", { + .structure <- function(pval_mat) { + id1 <- apply(pval_mat[, c(2, 4)], 1, which.min) + return(cbind( + pval_mat[, c(1, 3)][cbind(seq_along(id1), id1)], + pval_mat[, c(2, 4)][cbind(seq_along(id1), id1)]) + ) + } + for(p in c(0.01, 0.05, 0.1)) { + theta <- .Call("test_find_theta", test_pwm, snpInfo$prior, trans_mat, quantile(c(scores), 1 - p), package = "atSNP") + p_values <- .Call("test_p_value", test_pwm, snpInfo$prior, snpInfo$transition, c(scores), theta, 100, package = "atSNP") + p_values_s <- .structure(p_values) + expect_equal(p_values_s[, 2], apply(p_values[, c(2, 4)], 1, min)) + } + }) Test passed 😸 > > proc.time() user system elapsed 17.863 0.850 28.147
atSNP.Rcheck/atSNP-Ex.timings
name | user | system | elapsed | |
ComputeMotifScore | 0.991 | 0.123 | 1.623 | |
ComputePValues | 1.155 | 0.261 | 1.956 | |
GetIUPACSequence | 0.002 | 0.002 | 0.007 | |
LoadFastaData | 0.270 | 0.057 | 0.496 | |
LoadMotifLibrary | 0.261 | 0.035 | 0.452 | |
LoadSNPData | 0 | 0 | 0 | |
MatchSubsequence | 0.937 | 0.187 | 1.646 | |
dtMotifMatch | 0.982 | 0.260 | 1.814 | |
plotMotifMatch | 1.310 | 0.208 | 2.312 | |