diff --git a/R/ContingencyTableTest.R b/R/ContingencyTableTest.R
index 0029172e..0c9d9adc 100644
--- a/R/ContingencyTableTest.R
+++ b/R/ContingencyTableTest.R
@@ -36,7 +36,7 @@ ContingencyTableTest <- R6Class(
                 row_loc = rep.int(seq_len(r), row_sum) - 1,
                 col_loc = rep.int(seq_len(c), col_sum) - 1,
                 statistic_func = private$.statistic_func,
-                n_permu = as.integer(private$.n_permu)
+                n_permu = private$.n_permu
             )
         }
     )
diff --git a/R/KSampleTest.R b/R/KSampleTest.R
index 39e8b756..c4e073a4 100644
--- a/R/KSampleTest.R
+++ b/R/KSampleTest.R
@@ -40,7 +40,7 @@ KSampleTest <- R6Class(
                 data = unname(private$.data),
                 group = as.integer(names(private$.data)),
                 statistic_func = private$.statistic_func,
-                n_permu = as.integer(private$.n_permu)
+                n_permu = private$.n_permu
             )
         }
     )
diff --git a/R/MultipleComparison.R b/R/MultipleComparison.R
index ec852e9d..e080a014 100644
--- a/R/MultipleComparison.R
+++ b/R/MultipleComparison.R
@@ -52,7 +52,7 @@ MultipleComparison <- R6Class(
                 data = unname(private$.data),
                 group = as.integer(names(private$.data)),
                 statistic_func = private$.statistic_func,
-                n_permu = as.integer(private$.n_permu)
+                n_permu = private$.n_permu
             )
         },
 
diff --git a/R/RCBD.R b/R/RCBD.R
index 15f5fb44..3c1e15c3 100644
--- a/R/RCBD.R
+++ b/R/RCBD.R
@@ -36,7 +36,7 @@ RCBD <- R6Class(
             private$.statistic_permu <- rcbd_pmt(
                 data = apply(private$.data, 2, sort),
                 statistic_func = private$.statistic_func,
-                n_permu = as.integer(private$.n_permu)
+                n_permu = private$.n_permu
             )
         }
     )
diff --git a/R/TwoSampleAssociationTest.R b/R/TwoSampleAssociationTest.R
index e6b4faab..c7986df7 100644
--- a/R/TwoSampleAssociationTest.R
+++ b/R/TwoSampleAssociationTest.R
@@ -31,7 +31,7 @@ TwoSampleAssociationTest <- R6Class(
                 x = data_y_order$x,
                 y = data_y_order$y,
                 statistic_func = private$.statistic_func,
-                n_permu = as.integer(private$.n_permu)
+                n_permu = private$.n_permu
             )
         }
     )
diff --git a/R/TwoSamplePairedTest.R b/R/TwoSamplePairedTest.R
index 67be7f7b..c76c1db1 100644
--- a/R/TwoSamplePairedTest.R
+++ b/R/TwoSamplePairedTest.R
@@ -35,7 +35,7 @@ TwoSamplePairedTest <- R6Class(
             private$.statistic_permu <- paired_pmt(
                 n = nrow(private$.data),
                 statistic_func = private$.statistic_func,
-                n_permu = as.integer(private$.n_permu)
+                n_permu = private$.n_permu
             )
         }
     )
diff --git a/R/TwoSampleTest.R b/R/TwoSampleTest.R
index fb6f6181..8834883d 100644
--- a/R/TwoSampleTest.R
+++ b/R/TwoSampleTest.R
@@ -43,7 +43,7 @@ TwoSampleTest <- R6Class(
                 data = c(private$.data$x, private$.data$y),
                 where_y = rep.int(c(FALSE, TRUE), c(n_1, n_2)),
                 statistic_func = private$.statistic_func,
-                n_permu = as.integer(private$.n_permu)
+                n_permu = private$.n_permu
             )
         }
     )
diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp
index f0d3b5d8..eb7a3099 100644
--- a/src/RcppExports.cpp
+++ b/src/RcppExports.cpp
@@ -11,7 +11,7 @@ Rcpp::Rostream<false>& Rcpp::Rcerr = Rcpp::Rcpp_cerr_get();
 #endif
 
 // association_pmt
-NumericVector association_pmt(const NumericVector x, NumericVector y, const Function statistic_func, const unsigned n_permu);
+NumericVector association_pmt(const NumericVector x, NumericVector y, const Function statistic_func, const R_xlen_t n_permu);
 RcppExport SEXP _LearnNonparam_association_pmt(SEXP xSEXP, SEXP ySEXP, SEXP statistic_funcSEXP, SEXP n_permuSEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
@@ -19,13 +19,13 @@ BEGIN_RCPP
     Rcpp::traits::input_parameter< const NumericVector >::type x(xSEXP);
     Rcpp::traits::input_parameter< NumericVector >::type y(ySEXP);
     Rcpp::traits::input_parameter< const Function >::type statistic_func(statistic_funcSEXP);
-    Rcpp::traits::input_parameter< const unsigned >::type n_permu(n_permuSEXP);
+    Rcpp::traits::input_parameter< const R_xlen_t >::type n_permu(n_permuSEXP);
     rcpp_result_gen = Rcpp::wrap(association_pmt(x, y, statistic_func, n_permu));
     return rcpp_result_gen;
 END_RCPP
 }
 // ksample_pmt
-NumericVector ksample_pmt(const NumericVector data, IntegerVector group, const Function statistic_func, const unsigned n_permu);
+NumericVector ksample_pmt(const NumericVector data, IntegerVector group, const Function statistic_func, const R_xlen_t n_permu);
 RcppExport SEXP _LearnNonparam_ksample_pmt(SEXP dataSEXP, SEXP groupSEXP, SEXP statistic_funcSEXP, SEXP n_permuSEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
@@ -33,13 +33,13 @@ BEGIN_RCPP
     Rcpp::traits::input_parameter< const NumericVector >::type data(dataSEXP);
     Rcpp::traits::input_parameter< IntegerVector >::type group(groupSEXP);
     Rcpp::traits::input_parameter< const Function >::type statistic_func(statistic_funcSEXP);
-    Rcpp::traits::input_parameter< const unsigned >::type n_permu(n_permuSEXP);
+    Rcpp::traits::input_parameter< const R_xlen_t >::type n_permu(n_permuSEXP);
     rcpp_result_gen = Rcpp::wrap(ksample_pmt(data, group, statistic_func, n_permu));
     return rcpp_result_gen;
 END_RCPP
 }
 // multicomp_pmt
-NumericVector multicomp_pmt(const IntegerVector group_i, const IntegerVector group_j, const NumericVector data, IntegerVector group, const Function statistic_func, const unsigned n_permu);
+NumericVector multicomp_pmt(const IntegerVector group_i, const IntegerVector group_j, const NumericVector data, IntegerVector group, const Function statistic_func, const R_xlen_t n_permu);
 RcppExport SEXP _LearnNonparam_multicomp_pmt(SEXP group_iSEXP, SEXP group_jSEXP, SEXP dataSEXP, SEXP groupSEXP, SEXP statistic_funcSEXP, SEXP n_permuSEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
@@ -49,39 +49,39 @@ BEGIN_RCPP
     Rcpp::traits::input_parameter< const NumericVector >::type data(dataSEXP);
     Rcpp::traits::input_parameter< IntegerVector >::type group(groupSEXP);
     Rcpp::traits::input_parameter< const Function >::type statistic_func(statistic_funcSEXP);
-    Rcpp::traits::input_parameter< const unsigned >::type n_permu(n_permuSEXP);
+    Rcpp::traits::input_parameter< const R_xlen_t >::type n_permu(n_permuSEXP);
     rcpp_result_gen = Rcpp::wrap(multicomp_pmt(group_i, group_j, data, group, statistic_func, n_permu));
     return rcpp_result_gen;
 END_RCPP
 }
 // paired_pmt
-NumericVector paired_pmt(const unsigned n, const Function statistic_func, const unsigned n_permu);
+NumericVector paired_pmt(const R_len_t n, const Function statistic_func, const R_xlen_t n_permu);
 RcppExport SEXP _LearnNonparam_paired_pmt(SEXP nSEXP, SEXP statistic_funcSEXP, SEXP n_permuSEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
     Rcpp::RNGScope rcpp_rngScope_gen;
-    Rcpp::traits::input_parameter< const unsigned >::type n(nSEXP);
+    Rcpp::traits::input_parameter< const R_len_t >::type n(nSEXP);
     Rcpp::traits::input_parameter< const Function >::type statistic_func(statistic_funcSEXP);
-    Rcpp::traits::input_parameter< const unsigned >::type n_permu(n_permuSEXP);
+    Rcpp::traits::input_parameter< const R_xlen_t >::type n_permu(n_permuSEXP);
     rcpp_result_gen = Rcpp::wrap(paired_pmt(n, statistic_func, n_permu));
     return rcpp_result_gen;
 END_RCPP
 }
 // rcbd_pmt
-NumericVector rcbd_pmt(NumericMatrix data, const Function statistic_func, const unsigned n_permu);
+NumericVector rcbd_pmt(NumericMatrix data, const Function statistic_func, const R_xlen_t n_permu);
 RcppExport SEXP _LearnNonparam_rcbd_pmt(SEXP dataSEXP, SEXP statistic_funcSEXP, SEXP n_permuSEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
     Rcpp::RNGScope rcpp_rngScope_gen;
     Rcpp::traits::input_parameter< NumericMatrix >::type data(dataSEXP);
     Rcpp::traits::input_parameter< const Function >::type statistic_func(statistic_funcSEXP);
-    Rcpp::traits::input_parameter< const unsigned >::type n_permu(n_permuSEXP);
+    Rcpp::traits::input_parameter< const R_xlen_t >::type n_permu(n_permuSEXP);
     rcpp_result_gen = Rcpp::wrap(rcbd_pmt(data, statistic_func, n_permu));
     return rcpp_result_gen;
 END_RCPP
 }
 // table_pmt
-NumericVector table_pmt(IntegerVector row_loc, const IntegerVector col_loc, const Function statistic_func, const unsigned n_permu);
+NumericVector table_pmt(IntegerVector row_loc, const IntegerVector col_loc, const Function statistic_func, const R_xlen_t n_permu);
 RcppExport SEXP _LearnNonparam_table_pmt(SEXP row_locSEXP, SEXP col_locSEXP, SEXP statistic_funcSEXP, SEXP n_permuSEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
@@ -89,13 +89,13 @@ BEGIN_RCPP
     Rcpp::traits::input_parameter< IntegerVector >::type row_loc(row_locSEXP);
     Rcpp::traits::input_parameter< const IntegerVector >::type col_loc(col_locSEXP);
     Rcpp::traits::input_parameter< const Function >::type statistic_func(statistic_funcSEXP);
-    Rcpp::traits::input_parameter< const unsigned >::type n_permu(n_permuSEXP);
+    Rcpp::traits::input_parameter< const R_xlen_t >::type n_permu(n_permuSEXP);
     rcpp_result_gen = Rcpp::wrap(table_pmt(row_loc, col_loc, statistic_func, n_permu));
     return rcpp_result_gen;
 END_RCPP
 }
 // twosample_pmt
-NumericVector twosample_pmt(const NumericVector data, LogicalVector where_y, const Function statistic_func, const unsigned n_permu);
+NumericVector twosample_pmt(const NumericVector data, LogicalVector where_y, const Function statistic_func, const R_xlen_t n_permu);
 RcppExport SEXP _LearnNonparam_twosample_pmt(SEXP dataSEXP, SEXP where_ySEXP, SEXP statistic_funcSEXP, SEXP n_permuSEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
@@ -103,7 +103,7 @@ BEGIN_RCPP
     Rcpp::traits::input_parameter< const NumericVector >::type data(dataSEXP);
     Rcpp::traits::input_parameter< LogicalVector >::type where_y(where_ySEXP);
     Rcpp::traits::input_parameter< const Function >::type statistic_func(statistic_funcSEXP);
-    Rcpp::traits::input_parameter< const unsigned >::type n_permu(n_permuSEXP);
+    Rcpp::traits::input_parameter< const R_xlen_t >::type n_permu(n_permuSEXP);
     rcpp_result_gen = Rcpp::wrap(twosample_pmt(data, where_y, statistic_func, n_permu));
     return rcpp_result_gen;
 END_RCPP
diff --git a/src/association_pmt.cpp b/src/association_pmt.cpp
index eb2501ac..90df70ff 100644
--- a/src/association_pmt.cpp
+++ b/src/association_pmt.cpp
@@ -7,7 +7,7 @@ NumericVector association_pmt(
     const NumericVector x,
     NumericVector y,
     const Function statistic_func,
-    const unsigned n_permu)
+    const R_xlen_t n_permu)
 {
     auto association_statistic = [&]() -> double {
         return as<double>(statistic_func(x, y));
diff --git a/src/ksample_pmt.cpp b/src/ksample_pmt.cpp
index 3b39d662..75034f81 100644
--- a/src/ksample_pmt.cpp
+++ b/src/ksample_pmt.cpp
@@ -7,7 +7,7 @@ NumericVector ksample_pmt(
     const NumericVector data,
     IntegerVector group,
     const Function statistic_func,
-    const unsigned n_permu)
+    const R_xlen_t n_permu)
 {
     auto ksample_statistic = [&]() -> double {
         return as<double>(statistic_func(data, group));
diff --git a/src/multicomp_pmt.cpp b/src/multicomp_pmt.cpp
index c60e55a5..bd2c5c27 100644
--- a/src/multicomp_pmt.cpp
+++ b/src/multicomp_pmt.cpp
@@ -9,19 +9,19 @@ NumericVector multicomp_pmt(
     const NumericVector data,
     IntegerVector group,
     const Function statistic_func,
-    const unsigned n_permu)
+    const R_xlen_t n_permu)
 {
-    unsigned n_group = group[group.size() - 1];
-    unsigned n_pair = n_group * (n_group - 1) / 2;
+    R_len_t n_group = group[group.size() - 1];
+    R_len_t n_pair = n_group * (n_group - 1) / 2;
 
     List split(n_group);
     auto do_split = [&]() {
-        for (unsigned i = 1; i <= n_group; i++) {
+        for (R_len_t i = 1; i <= n_group; i++) {
             split[i - 1] = data[group == i];
         }
     };
 
-    unsigned j;
+    R_len_t j;
 
     auto multicomp_statistic = [&]() -> double {
         return as<double>(statistic_func(split[group_i[j]], split[group_j[j]], data, group));
diff --git a/src/paired_pmt.cpp b/src/paired_pmt.cpp
index caa23619..1d85d492 100644
--- a/src/paired_pmt.cpp
+++ b/src/paired_pmt.cpp
@@ -4,24 +4,23 @@ using namespace Rcpp;
 
 // [[Rcpp::export]]
 NumericVector paired_pmt(
-    const unsigned n,
+    const R_len_t n,
     const Function statistic_func,
-    const unsigned n_permu)
+    const R_xlen_t n_permu)
 {
     LogicalVector swapped(n);
 
-    unsigned i = 0;
+    R_xlen_t i = 0;
+    R_xlen_t total = (1 << n);
 
     auto paired_statistic = [&]() -> double {
-        for (unsigned j = 0; j < n; j++) {
+        for (R_len_t j = 0; j < n; j++) {
             swapped[j] = ((i & (1 << j)) != 0);
         }
 
         return as<double>(statistic_func(swapped));
     };
 
-    unsigned total = (1 << n);
-
     if (n_permu == 0) {
         PermuBar bar(total, true);
 
diff --git a/src/rcbd_pmt.cpp b/src/rcbd_pmt.cpp
index b73c524c..1803930e 100644
--- a/src/rcbd_pmt.cpp
+++ b/src/rcbd_pmt.cpp
@@ -6,17 +6,17 @@ using namespace Rcpp;
 NumericVector rcbd_pmt(
     NumericMatrix data,
     const Function statistic_func,
-    const unsigned n_permu)
+    const R_xlen_t n_permu)
 {
     auto rcbd_statistic = [&]() -> double {
         return as<double>(statistic_func(data));
     };
 
-    unsigned i = 0;
-    unsigned n_col = data.ncol();
+    R_len_t i = 0;
+    R_len_t n_col = data.ncol();
     if (n_permu == 0) {
-        unsigned total = 1;
-        for (unsigned j = 0; j < n_col; j++) {
+        R_xlen_t total = 1;
+        for (R_len_t j = 0; j < n_col; j++) {
             total *= n_permutation(data.column(j));
         }
 
diff --git a/src/table_pmt.cpp b/src/table_pmt.cpp
index c2ba78c3..a6e59c71 100644
--- a/src/table_pmt.cpp
+++ b/src/table_pmt.cpp
@@ -7,15 +7,15 @@ NumericVector table_pmt(
     IntegerVector row_loc,
     const IntegerVector col_loc,
     const Function statistic_func,
-    const unsigned n_permu)
+    const R_xlen_t n_permu)
 {
-    unsigned n = row_loc.size();
+    R_len_t n = row_loc.size();
 
     IntegerMatrix data(no_init(row_loc[n - 1] + 1, col_loc[n - 1] + 1));
 
     auto table_statistic = [&]() -> double {
         data.fill(0);
-        for (unsigned i = 0; i < n; i++) {
+        for (R_len_t i = 0; i < n; i++) {
             data(row_loc[i], col_loc[i])++;
         }
 
diff --git a/src/twosample_pmt.cpp b/src/twosample_pmt.cpp
index 6dfc637a..55e34d6b 100644
--- a/src/twosample_pmt.cpp
+++ b/src/twosample_pmt.cpp
@@ -7,7 +7,7 @@ NumericVector twosample_pmt(
     const NumericVector data,
     LogicalVector where_y,
     const Function statistic_func,
-    const unsigned n_permu)
+    const R_xlen_t n_permu)
 {
     auto twosample_statistic = [&]() -> double {
         return as<double>(statistic_func(data[!where_y], data[where_y]));
diff --git a/src/utils.h b/src/utils.h
index a0197586..6cce7df2 100644
--- a/src/utils.h
+++ b/src/utils.h
@@ -10,7 +10,8 @@ using namespace Rcpp;
 
 // runif(max = n) (tied to the same RNG which R uses)
 
-inline unsigned rand_int(const unsigned& n)
+template <typename T>
+T rand_int(T n)
 {
     return floor(unif_rand() * n);
 }
@@ -20,9 +21,9 @@ inline unsigned rand_int(const unsigned& n)
 template <typename T>
 void random_shuffle(T&& v)
 {
-    unsigned j;
-    unsigned n = v.size();
-    for (unsigned i = 0; i < n - 1; i++) {
+    R_len_t j;
+    R_len_t n = v.size();
+    for (R_len_t i = 0; i < n - 1; i++) {
         j = i + rand_int(n - i);
         std::swap(v[i], v[j]);
     }
@@ -31,14 +32,14 @@ void random_shuffle(T&& v)
 // number of permutations
 
 template <typename T>
-unsigned n_permutation(T&& v)
+R_xlen_t n_permutation(T&& v)
 {
     double A = 1;
 
-    unsigned n_i = 0;
-    unsigned n = v.size();
+    R_len_t n_i = 0;
+    R_len_t n = v.size();
     double current = v[0];
-    for (unsigned i = 0; i < n; i++) {
+    for (R_len_t i = 0; i < n; i++) {
         A *= (i + 1);
         if (v[i] == current) {
             n_i++;
@@ -49,17 +50,17 @@ unsigned n_permutation(T&& v)
         current = v[i];
     }
 
-    return (unsigned)A;
+    return (R_xlen_t)A;
 }
 
 // progress bar
 
 class PermuBar {
 private:
-    unsigned _total;
+    R_xlen_t _total;
 
-    unsigned _update_i = 0;
-    unsigned _update_every;
+    R_xlen_t _update_i = 0;
+    R_xlen_t _update_every;
 
     NumericVector::iterator _iter;
     NumericVector::iterator _end;
@@ -74,7 +75,7 @@ class PermuBar {
 
         buffer << "\015";
 
-        unsigned percent = 100 - 100 * (_end - _iter) / _total;
+        unsigned percent = 100 - (unsigned)(100 * (_end - _iter) / _total);
 
         buffer << "\033[31m" << percent << "%";
 
@@ -103,7 +104,7 @@ class PermuBar {
 public:
     NumericVector statistic_permu;
 
-    PermuBar(unsigned n_permu, bool exact, unsigned statistic_size = 1)
+    PermuBar(R_xlen_t n_permu, bool exact, R_len_t statistic_size = 1)
     {
         _total = n_permu * statistic_size;