Merge pull request #178 from kgoldfeld/fix-addCorGen-2

Fix add cor gen 2

Author: kgoldfeld

CITATION.cff

@@ -123,25 +123,6 @@ references:
     email: matteo.fasiolo@gmail.com
   year: '2022'
   url: https://CRAN.R-project.org/package=mvnfast
-- type: software
-  title: mvtnorm
-  abstract: 'mvtnorm: Multivariate Normal and t Distributions'
-  notes: Imports
-  authors:
-  - family-names: Genz
-    given-names: Alan
-  - family-names: Bretz
-    given-names: Frank
-  - family-names: Miwa
-    given-names: Tetsuhisa
-  - family-names: Mi
-    given-names: Xuefei
-  - family-names: Hothorn
-    given-names: Torsten
-    email: Torsten.Hothorn@R-project.org
-    orcid: https://orcid.org/0000-0001-8301-0471
-  year: '2022'
-  url: https://CRAN.R-project.org/package=mvtnorm
 - type: software
   title: Rcpp
   abstract: 'Rcpp: Seamless R and C++ Integration'
@@ -498,6 +479,16 @@ references:
     given-names: Emanuele
   year: '2022'
   url: https://CRAN.R-project.org/package=blockmatrix
+- type: software
+  title: pbv
+  abstract: 'pbv: Probabilities for Bivariate Normal Distribution'
+  notes: LinkingTo
+  authors:
+  - family-names: Robitzsch
+    given-names: Alexander
+  year: '2022'
+  url: https://CRAN.R-project.org/package=pbv
+  version: '>= 0.4-22'
 identifiers:
 - type: url
   value: https://kgoldfeld.github.io/simstudy/dev/

DESCRIPTION

@@ -2,7 +2,7 @@ Type: Package
 Package: simstudy
 Title: Simulation of Study Data
 Version: 0.5.1.9000
-Date: 2022-10-03
+Date: 2022-11-16
 Authors@R: 
     c(person(given = "Keith",
              family = "Goldfeld",
@@ -33,7 +33,6 @@ Imports:
     glue,
     methods,
     mvnfast,
-    mvtnorm,
     Rcpp,
     backports
 Suggests: 
@@ -60,8 +59,9 @@ Suggests:
     survminer,
     blockmatrix
 LinkingTo: 
-    Rcpp
+    Rcpp, 
+    pbv (>= 0.4-22)
 VignetteBuilder: 
     knitr
 Encoding: UTF-8
-RoxygenNote: 7.2.1
+RoxygenNote: 7.2.2

NEWS.md

@@ -5,6 +5,12 @@
 accommodate clustered observations over time where the within-cluster correlation
 in the same time period can be different from the within-cluster correlation across time periods.
 
+## Major fixes
+
+* Overhauled function `addCorGen` to make it more flexible. It can now handle cluster
+dependent data, and not just time-dependent data. In addition, performance has been
+dramatically improved.
+
 # simstudy 0.5.1
 
 ## Minor fixes

R/RcppExports.R

@@ -17,3 +17,15 @@ chkNonIncreasing <- function(adjmatrix) {
     .Call(`_simstudy_chkNonIncreasing`, adjmatrix)
 }
 
+checkBoundsBin <- function(p1, p2, d) {
+    invisible(.Call(`_simstudy_checkBoundsBin`, p1, p2, d))
+}
+
+findRhoBin <- function(p1, p2, d) {
+    .Call(`_simstudy_findRhoBin`, p1, p2, d)
+}
+
+getRhoMat <- function(N, P, TCORR) {
+    .Call(`_simstudy_getRhoMat`, N, P, TCORR)
+}
+

R/add_correlated_data.R

@@ -258,36 +258,21 @@ genCorFlex <- function(n, defs, rho = 0, tau = NULL, corstr = "cs", corMatrix =
 #' Multivariate Binary Variates. The American Statistician 1991;45:302-4.
 #' @examples
 #' set.seed(23432)
-#' l <- c(8, 10, 12)
+#' lambda <- c(8, 10, 12)
 #'
-#' genCorGen(1000, nvars = 3, params1 = l, dist = "poisson", rho = .7, corstr = "cs")
-#' genCorGen(1000, nvars = 3, params1 = 5, dist = "poisson", rho = .7, corstr = "cs")
-#' genCorGen(1000, nvars = 3, params1 = l, dist = "poisson", rho = .7, corstr = "cs", wide = TRUE)
-#' genCorGen(1000, nvars = 3, params1 = 5, dist = "poisson", rho = .7, corstr = "cs", wide = TRUE)
+#' genCorGen(100, nvars = 3, params1 = lambda, dist = "poisson", rho = .7, corstr = "cs")
+#' genCorGen(100, nvars = 3, params1 = 5, dist = "poisson", rho = .7, corstr = "cs")
+#' genCorGen(100, nvars = 3, params1 = lambda, dist = "poisson", rho = .7, corstr = "cs", wide = TRUE)
+#' genCorGen(100, nvars = 3, params1 = 5, dist = "poisson", rho = .7, corstr = "cs", wide = TRUE)
 #'
-#' genCorGen(1000,
-#'   nvars = 3, params1 = l, dist = "poisson", rho = .7, corstr = "cs",
+#' genCorGen(100,
+#'   nvars = 3, params1 = lambda, dist = "poisson", rho = .7, corstr = "cs",
 #'   cnames = "new_var"
 #' )
-#' genCorGen(1000,
-#'   nvars = 3, params1 = l, dist = "poisson", rho = .7, corstr = "cs",
+#' genCorGen(100,
+#'   nvars = 3, params1 = lambda, dist = "poisson", rho = .7, corstr = "cs",
 #'   wide = TRUE, cnames = "a, b, c"
 #' )
-#'
-#' genCorGen(1000, nvars = 3, params1 = c(.3, .5, .7), dist = "binary", rho = .3, corstr = "cs")
-#' genCorGen(1000,
-#'   nvars = 3, params1 = l, params2 = c(1, 1, 1), dist = "gamma", rho = .3,
-#'   corstr = "cs", wide = TRUE
-#' )
-#'
-#' genCorGen(1000,
-#'   nvars = 3, params1 = c(.3, .5, .7), dist = "binary",
-#'   corMatrix = genCorMat(3), method = "ep"
-#' )
-#' genCorGen(1000,
-#'   nvars = 3, params1 = c(.3, .5, .7), dist = "binary",
-#'   corMatrix = genCorMat(3), method = "copula"
-#' )
 #' @export
 #' @concept correlated
 genCorGen <- function(n, nvars, params1, params2 = NULL, dist, rho, corstr,
@@ -359,6 +344,10 @@ genCorGen <- function(n, nvars, params1, params2 = NULL, dist, rho, corstr,
   if (dist != "binary" & method == "ep") {
     stop("Method `ep` applies only to binary data generation")
   }
+  
+  if (!is.null(corMatrix)) {
+    assertClass(corMatrix = corMatrix, class = "matrix")
+  }
 
   ####
 
@@ -404,7 +393,7 @@ genCorGen <- function(n, nvars, params1, params2 = NULL, dist, rho, corstr,
 
     if (!is.null(cnames)) setnames(dFinal, "X", cnames)
   } else {
-    dFinal <- dcast(dtM, id ~ seq, value.var = "X")
+    dFinal <- data.table::dcast(dtM, id ~ seq, value.var = "X")
     if (!is.null(cnames)) {
       nnames <- trimws(unlist(strsplit(cnames, split = ",")))
       setnames(dFinal, paste0("V", 1:nvars), nnames)
@@ -416,66 +405,6 @@ genCorGen <- function(n, nvars, params1, params2 = NULL, dist, rho, corstr,
   return(dFinal[])
 }
 
-
-# TODO Implement Emrich and Piedmonte algorithm for correlated binary data?
-#' Internal functions called by genCorGen and addCorGen - returns matrix
-#'
-#' @param nvars Number of new variables to generate
-#' @param corMatrix Correlation matrix
-#' @param rho Correlation coefficient
-#' @param corstr Correlation structure
-#' @return A correlation matrix
-#' @noRd
-.checkBoundsBin <- function(p1, p2, d) {
-  l <- (p1 * p2) / ((1 - p1) * (1 - p2))
-  L <- max(-sqrt(l), -sqrt(1 / l))
-
-  u <- (p1 * (1 - p2)) / (p2 * (1 - p1))
-  U <- min(sqrt(u), sqrt(1 / u))
-
-  if ((d < L & isTRUE(all.equal(d, L)) == FALSE) |
-    (d > U & isTRUE(all.equal(d, U)) == FALSE)) {
-    LU <- paste0("(", round(L, 2), " ... ", round(U, 2), ")")
-    stopText <- paste("Specified correlation", d, "out of range", LU)
-    stop(stopText)
-  }
-}
-
-
-#'
-.findRhoBin <- function(p1, p2, d) {
-  .checkBoundsBin(p1, p2, d)
-
-  target <- d * sqrt(p1 * p2 * (1 - p1) * (1 - p2)) + p1 * p2
-
-  # given p1, p2 & d, bisection search for corresponding rho
-
-  Max <- 1
-  Min <- -1
-  test <- 0
-  found <- FALSE
-
-  while (!found) {
-    corr <- diag(2)
-    corr[1, 2] <- corr[2, 1] <- test
-
-    est <- mvtnorm::pmvnorm(lower = rep(-Inf, 2), upper = c(stats::qnorm(p1), stats::qnorm(p2)), mean = c(0, 0), corr = corr)
-
-    if (round(est, 5) == round(target, 5)) {
-      found <- TRUE
-      rho <- test
-    } else if (est < target) {
-      Min <- test
-      test <- (Min + Max) / 2
-    } else {
-      Max <- test
-      test <- (Min + Max) / 2
-    }
-  }
-
-  return(rho)
-}
-
 #'
 .genBinEP <- function(n, p, tcorr) {
 
@@ -486,16 +415,9 @@ genCorGen <- function(n, nvars, params1, params2 = NULL, dist, rho, corstr,
 
   np <- length(p)
   phicorr <- diag(length(p))
-
-  for (i in (1:(np - 1))) {
-    for (j in ((i + 1):np)) {
-      p1 <- p[i]
-      p2 <- p[j]
-
-      phicorr[j, i] <- phicorr[i, j] <- .findRhoBin(p1, p2, tcorr[i, j])
-    }
-  }
-
+  
+  phicorr <- getRhoMat(np, p, tcorr)
+  
   # check that phicorr is positive definite (PD), if not adjust to nearest PD matrix
   if (!all(eigen(phicorr)$values > 0)) {
     phicorr <- Matrix::nearPD(phicorr)$mat
@@ -632,6 +554,8 @@ genBlockMat <- function(rho, nInds, nPeriods, corstr = "ind",
   assertAtLeast(nPeriods = nPeriods, minVal = 2)
   assertInRange(rho = rho, range = c(-1,1))
 
+  assertOption(corstr = corstr, options = c("ind", "cs", "ar1"))
+  
   if (!is.null(iRho)) {
     assertInRange(iRho = iRho, range = c(-1,1))
     if (length(rho) == 1) {
@@ -647,7 +571,6 @@ genBlockMat <- function(rho, nInds, nPeriods, corstr = "ind",
     }
   }
 
-  assertOption(corstr = corstr, options = c("ind", "cs", "ar1"))
 
   ###
 

R/generate_correlated_data.R

@@ -652,7 +652,7 @@ viewBasis <- function(knots, degree) {
   )
 
   ggplot2::ggplot(data = dtmelt, ggplot2::aes(x = x, y = value, group = basis)) +
-    ggplot2::geom_line(ggplot2::aes(color = basis), size = 1) +
+    ggplot2::geom_line(ggplot2::aes(color = basis), linewidth = 1) +
     ggplot2::theme(
       legend.position = "none",
       panel.grid.minor = ggplot2::element_blank()
@@ -727,7 +727,8 @@ viewSplines <- function(knots, degree, theta) {
   dx[, Spline := factor(index)]
 
   p <- ggplot2::ggplot(data = dx) +
-    ggplot2::geom_line(ggplot2::aes(x = x, y = y.spline, color = Spline), size = 1) +
+    ggplot2::geom_line(ggplot2::aes(x = x, y = y.spline, color = Spline), 
+        linewidth = 1) +
     ggplot2::scale_y_continuous(limits = c(0, 1)) +
     ggplot2::scale_x_continuous(limits = c(0, 1), breaks = knots) +
     ggplot2::theme(panel.grid.minor = ggplot2::element_blank()) +
@@ -841,7 +842,7 @@ survParamPlot <- function(formula, shape, points = NULL, n = 100, scale = 1) {
   )
 
   p <- ggplot2::ggplot(data = dd, ggplot2::aes(x = T, y = p)) +
-    ggplot2::geom_line(size = 0.8) +
+    ggplot2::geom_line(linewidth = 0.8) +
     ggplot2::scale_y_continuous(limits = c(0,1), name = "probability of survival") +
     ggplot2::scale_x_continuous(name = "time") +
     ggplot2::theme(panel.grid = ggplot2::element_blank(),

R/utility.R

@@ -11,18 +11,18 @@ Rcpp::Rostream<false>& Rcpp::Rcerr = Rcpp::Rcpp_cerr_get();
 #endif
 
 // matMultinom
-Rcpp::IntegerVector matMultinom(Rcpp::NumericMatrix probmatrix);
+IntegerVector matMultinom(NumericMatrix probmatrix);
 RcppExport SEXP _simstudy_matMultinom(SEXP probmatrixSEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
     Rcpp::RNGScope rcpp_rngScope_gen;
-    Rcpp::traits::input_parameter< Rcpp::NumericMatrix >::type probmatrix(probmatrixSEXP);
+    Rcpp::traits::input_parameter< NumericMatrix >::type probmatrix(probmatrixSEXP);
     rcpp_result_gen = Rcpp::wrap(matMultinom(probmatrix));
     return rcpp_result_gen;
 END_RCPP
 }
 // markovChains
-Rcpp::IntegerMatrix markovChains(int nchains, NumericMatrix P, int chainLen, IntegerVector state0);
+IntegerMatrix markovChains(int nchains, NumericMatrix P, int chainLen, IntegerVector state0);
 RcppExport SEXP _simstudy_markovChains(SEXP nchainsSEXP, SEXP PSEXP, SEXP chainLenSEXP, SEXP state0SEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
@@ -36,7 +36,7 @@ BEGIN_RCPP
 END_RCPP
 }
 // clipVec
-Rcpp::IntegerVector clipVec(IntegerVector id, IntegerVector seq, IntegerVector event);
+IntegerVector clipVec(IntegerVector id, IntegerVector seq, IntegerVector event);
 RcppExport SEXP _simstudy_clipVec(SEXP idSEXP, SEXP seqSEXP, SEXP eventSEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
@@ -49,22 +49,63 @@ BEGIN_RCPP
 END_RCPP
 }
 // chkNonIncreasing
-bool chkNonIncreasing(Rcpp::NumericMatrix adjmatrix);
+bool chkNonIncreasing(NumericMatrix adjmatrix);
 RcppExport SEXP _simstudy_chkNonIncreasing(SEXP adjmatrixSEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
     Rcpp::RNGScope rcpp_rngScope_gen;
-    Rcpp::traits::input_parameter< Rcpp::NumericMatrix >::type adjmatrix(adjmatrixSEXP);
+    Rcpp::traits::input_parameter< NumericMatrix >::type adjmatrix(adjmatrixSEXP);
     rcpp_result_gen = Rcpp::wrap(chkNonIncreasing(adjmatrix));
     return rcpp_result_gen;
 END_RCPP
 }
+// checkBoundsBin
+void checkBoundsBin(double p1, double p2, double d);
+RcppExport SEXP _simstudy_checkBoundsBin(SEXP p1SEXP, SEXP p2SEXP, SEXP dSEXP) {
+BEGIN_RCPP
+    Rcpp::RNGScope rcpp_rngScope_gen;
+    Rcpp::traits::input_parameter< double >::type p1(p1SEXP);
+    Rcpp::traits::input_parameter< double >::type p2(p2SEXP);
+    Rcpp::traits::input_parameter< double >::type d(dSEXP);
+    checkBoundsBin(p1, p2, d);
+    return R_NilValue;
+END_RCPP
+}
+// findRhoBin
+double findRhoBin(double p1, double p2, double d);
+RcppExport SEXP _simstudy_findRhoBin(SEXP p1SEXP, SEXP p2SEXP, SEXP dSEXP) {
+BEGIN_RCPP
+    Rcpp::RObject rcpp_result_gen;
+    Rcpp::RNGScope rcpp_rngScope_gen;
+    Rcpp::traits::input_parameter< double >::type p1(p1SEXP);
+    Rcpp::traits::input_parameter< double >::type p2(p2SEXP);
+    Rcpp::traits::input_parameter< double >::type d(dSEXP);
+    rcpp_result_gen = Rcpp::wrap(findRhoBin(p1, p2, d));
+    return rcpp_result_gen;
+END_RCPP
+}
+// getRhoMat
+NumericMatrix getRhoMat(int N, NumericVector P, NumericMatrix TCORR);
+RcppExport SEXP _simstudy_getRhoMat(SEXP NSEXP, SEXP PSEXP, SEXP TCORRSEXP) {
+BEGIN_RCPP
+    Rcpp::RObject rcpp_result_gen;
+    Rcpp::RNGScope rcpp_rngScope_gen;
+    Rcpp::traits::input_parameter< int >::type N(NSEXP);
+    Rcpp::traits::input_parameter< NumericVector >::type P(PSEXP);
+    Rcpp::traits::input_parameter< NumericMatrix >::type TCORR(TCORRSEXP);
+    rcpp_result_gen = Rcpp::wrap(getRhoMat(N, P, TCORR));
+    return rcpp_result_gen;
+END_RCPP
+}
 
 static const R_CallMethodDef CallEntries[] = {
     {"_simstudy_matMultinom", (DL_FUNC) &_simstudy_matMultinom, 1},
     {"_simstudy_markovChains", (DL_FUNC) &_simstudy_markovChains, 4},
     {"_simstudy_clipVec", (DL_FUNC) &_simstudy_clipVec, 3},
     {"_simstudy_chkNonIncreasing", (DL_FUNC) &_simstudy_chkNonIncreasing, 1},
+    {"_simstudy_checkBoundsBin", (DL_FUNC) &_simstudy_checkBoundsBin, 3},
+    {"_simstudy_findRhoBin", (DL_FUNC) &_simstudy_findRhoBin, 3},
+    {"_simstudy_getRhoMat", (DL_FUNC) &_simstudy_getRhoMat, 3},
     {NULL, NULL, 0}
 };
 

man/addCorGen.Rd

@@ -1,6 +1,11 @@
 #include <Rcpp.h>
-using namespace Rcpp;
 
+// [[Rcpp::depends(pbv)]]
+#include <pbv.h>
+
+#include<string.h>
+using namespace std;
+using namespace Rcpp;
 
 int vecMultinom(NumericVector probs) {
 
@@ -19,11 +24,11 @@ int vecMultinom(NumericVector probs) {
 }
 
 // [[Rcpp::export]]
-Rcpp::IntegerVector matMultinom(Rcpp::NumericMatrix probmatrix) {
+IntegerVector matMultinom(NumericMatrix probmatrix) {
 
   int rows = probmatrix.nrow();
 
-  Rcpp::IntegerVector ans(rows);
+  IntegerVector ans(rows);
 
   for(int i = 0; i < rows; ++i) {
     ans[i] = vecMultinom(probmatrix(i, _));
@@ -32,7 +37,7 @@ Rcpp::IntegerVector matMultinom(Rcpp::NumericMatrix probmatrix) {
   return(ans);
 }
 
-Rcpp::IntegerVector singleMarkovChain(NumericMatrix P, 
+IntegerVector singleMarkovChain(NumericMatrix P, 
                                       int chainLen, int state0) {
 
   IntegerVector states(chainLen);
@@ -46,7 +51,7 @@ Rcpp::IntegerVector singleMarkovChain(NumericMatrix P,
 }
 
 // [[Rcpp::export]]
-Rcpp::IntegerMatrix markovChains(int nchains, NumericMatrix P, 
+IntegerMatrix markovChains(int nchains, NumericMatrix P, 
                                  int chainLen, IntegerVector state0) {
 
   IntegerMatrix stateMat(nchains, chainLen);
@@ -60,7 +65,8 @@ Rcpp::IntegerMatrix markovChains(int nchains, NumericMatrix P,
 }
 
 // [[Rcpp::export]]
-Rcpp::IntegerVector clipVec(IntegerVector id, IntegerVector seq, IntegerVector event) {
+IntegerVector clipVec(IntegerVector id, IntegerVector seq, 
+                            IntegerVector event) {
 
   int uid = unique(id).length();
   IntegerVector last(uid);
@@ -87,7 +93,7 @@ Rcpp::IntegerVector clipVec(IntegerVector id, IntegerVector seq, IntegerVector e
 }
 
 // [[Rcpp::export]]
-bool chkNonIncreasing(Rcpp::NumericMatrix adjmatrix) {
+bool chkNonIncreasing(NumericMatrix adjmatrix) {
 
   int nr = adjmatrix.nrow(), nc = adjmatrix.ncol();
   bool unsorted = FALSE;
@@ -109,4 +115,94 @@ bool chkNonIncreasing(Rcpp::NumericMatrix adjmatrix) {
   return(unsorted);
 }
 
+// [[Rcpp::export]]
+void checkBoundsBin(double p1, double p2, double d) {
+  
+  double l =  (p1 * p2) / ((1 - p1) * (1 - p2));
+  double L = std::max(-sqrt(l), -sqrt(1 / l));
+  
+  double u = (p1 * (1 - p2)) / (p2 * (1 - p1));
+  double U = std::min(sqrt(u), sqrt(1 / u));
+  
+  string LU;
+  string stopText;
+  
+  if ((d < L)  | (d > U)) {
+    LU = "(" + to_string(L) + " ... " + to_string(U) + ")";
+    stopText = "Specified correlation " + to_string(d) + " out of range " + LU;
+    stop(stopText);
+  } 
+}
+
+// [[Rcpp::export]]
+double findRhoBin(double p1, double p2, double d) {
+  
+  checkBoundsBin(p1, p2, d);
+
+  double target;
+  double Max = 1;
+  double Min = -1;
+  double test = 0;
+  bool found = FALSE;
+  NumericVector bounds(2);
+  double est;
+  double rho;
+  NumericVector check(2);
+
+ target = d * sqrt(p1 * p2 * (1 - p1) * (1 - p2)) + p1 * p2;
+ bounds(0) = R::qnorm(p1, 0, 1, TRUE, FALSE);
+ bounds(1) = R::qnorm(p2, 0, 1, TRUE, FALSE); 
+  
+// given p1, p2 & d, bisection search for corresponding rho
+    
+  while (!found) {
+
+      est = pbv::pbv_rcpp_pbvnorm0(bounds(0), bounds(1), test);
+      
+      check(0) = est;
+      check(1) = target;
+      check = round(check, 5);
+      
+      if (check(0) == check(1)) {
+        found = TRUE;
+        rho = test;
+      } 
+      else if (est < target) {
+        Min = test;
+        test = (Min + Max) / 2;
+      } 
+      else {
+        Max = test;
+        test = (Min + Max) / 2;
+      } 
+    }
+    
+    return(rho);
+} 
+
+// [[Rcpp::export]]
+NumericMatrix getRhoMat(int N, NumericVector P, NumericMatrix TCORR) {
+  
+  NumericMatrix PCORR(TCORR.nrow(), TCORR.ncol());
+  double p1;
+  double p2;
+  double rho;
+  
+  for (int i = 0; i < (N - 1); i++) {
+    for (int j = (i+1); j < N; j++) {
+      p1 = P[i];
+      p2 = P[j];
+      rho = findRhoBin(p1, p2, TCORR(i, j));
+      PCORR(i, j) = rho;
+      PCORR(j, i) = rho;
+    }
+  }
+  
+  for (int i = 0; i < N; i++) {
+    PCORR(i, i) = 1;
+  }
+  
+  return(PCORR);
+  
+}
 

man/genCorGen.Rd

@@ -1,16 +1,15 @@
 # .checkBoundsBin ----
 test_that("Correlation boundaries for binary variables are correct", {
-  x1 <- c(0, 0, 0, 1, 0)
-  x2 <- c(1, 1, 1, 0, 0)
-
-  p1 <- mean(x1)
-  p2 <- mean(x2)
+  
+  p1 <- .5
+  p2 <- .8
 
-  expect_error(.checkBoundsBin(p1, p2, -1))
+  expect_error(checkBoundsBin(p1, p2, d = .9))
+  expect_error(checkBoundsBin(p1, p2, d = -.6))
 
-  expect_silent(.checkBoundsBin(p1, p2, -0.6))
-  expect_silent(.checkBoundsBin(p1, p2, 0.4))
-  expect_silent(.checkBoundsBin(p1, p2, cor(x1, x2)))
+  expect_silent(checkBoundsBin(p1, p2, -0.4))
+  expect_silent(checkBoundsBin(p1, p2, 0.3))
+  expect_silent(checkBoundsBin(p1, p2, 0.2))
 })
 
 # .findRhoBin ----