diff --git a/src/IC_SSH.cpp b/src/IC_SSH.cpp index cc24e5a..a744b9b 100644 --- a/src/IC_SSH.cpp +++ b/src/IC_SSH.cpp @@ -68,23 +68,26 @@ std::vector IC_SSHICM(const std::vector& d, // Step 2: Generate random permutations of s and compute IC for each std::vector IC_results(permutation_number, 0.0); // Store IC values for each permutation - // Initialize the random number generator with the seed - std::mt19937 gen(seed); + // Step 3: Generate a random seed using the input seed + std::mt19937 seed_gen(seed); // Initialize random number generator with the input seed + // std::uniform_int_distribution<> dis(1, std::numeric_limits::max()); // Define a distribution for random integers + std::uniform_int_distribution<> dis(1, 100); + int randomseed = dis(seed_gen); // Generate a random integer using the input seed - // Step 3: Perform parallel computation + // Step 4: Perform parallel computation RcppThread::parallelFor(0, permutation_number, [&](size_t i) { - // Generate a unique seed for each permutation - std::mt19937 local_gen(seed + i); // Modify seed for each thread + // Step 4.1: Generate a unique seed for each permutation by adding the iteration index to the randomseed + std::mt19937 local_gen(randomseed + i); // Modify seed for each thread - // Step 3.1: Permute d + // Step 4.2: Permute d std::vector permuted_d = d; // Copy the original d std::shuffle(permuted_d.begin(), permuted_d.end(), local_gen); // Shuffle based on the unique seed for each thread - // Step 3.2: Compute IC for the permuted d + // Step 4.3: Compute IC for the permuted d IC_results[i] = IC_SSH(permuted_d, s, bin_method); }); - // Step 4: Compute p-value by comparing permuted IC values to the true IC value + // Step 5: Compute p-value by comparing permuted IC values to the true IC value int greater_count = 0; for (size_t i = 0; i < IC_results.size(); ++i) { if (IC_results[i] >= true_IC) { @@ -98,36 +101,40 @@ std::vector IC_SSHICM(const std::vector& d, return {true_IC, p_value}; } - // // IC_SSHICM: Parallel computation of IC_SSH over permutations, returning IC value and p-value // // [[Rcpp::export]] // std::vector IC_SSHICM(const std::vector& d, // const std::vector& s, -// const std::vector>& permutation, -// const std::string& bin_method = "SquareRoot") { +// unsigned int seed, +// int permutation_number, +// const std::string& bin_method = "Sturges") { // if (s.size() != d.size()) { // throw std::invalid_argument("Vectors s and d must have the same length."); // } // -// size_t num_permutations = permutation.size(); // Number of permutations -// std::vector IC_results(num_permutations, 0.0); // Store IC values for each permutation -// -// // Calculate the true IC value using the original d and s +// // Step 1: Calculate the true IC value using the original d and s // double true_IC = IC_SSH(d, s, bin_method); // -// // Parallel computation using RcppThread -// RcppThread::parallelFor(0, num_permutations, [&](size_t i) { -// // Extract a single permutation from the matrix -// std::vector permuted_s(s.size()); -// for (size_t j = 0; j < s.size(); ++j) { -// permuted_s[j] = permutation[i][j]; // Permute s based on the current permutation -// } +// // Step 2: Generate random permutations of s and compute IC for each +// std::vector IC_results(permutation_number, 0.0); // Store IC values for each permutation +// +// // Initialize the random number generator with the seed +// std::mt19937 gen(seed); +// +// // Step 3: Perform parallel computation +// RcppThread::parallelFor(0, permutation_number, [&](size_t i) { +// // Generate a unique seed for each permutation +// std::mt19937 local_gen(seed + i); // Modify seed for each thread +// +// // Step 3.1: Permute d +// std::vector permuted_d = d; // Copy the original d +// std::shuffle(permuted_d.begin(), permuted_d.end(), local_gen); // Shuffle based on the unique seed for each thread // -// // Compute IC for the current permutation -// IC_results[i] = IC_SSH(d, permuted_s, bin_method); +// // Step 3.2: Compute IC for the permuted d +// IC_results[i] = IC_SSH(permuted_d, s, bin_method); // }); // -// // Compute p-value by comparing permuted IC values to the true IC value +// // Step 4: Compute p-value by comparing permuted IC values to the true IC value // int greater_count = 0; // for (size_t i = 0; i < IC_results.size(); ++i) { // if (IC_results[i] >= true_IC) { @@ -135,7 +142,7 @@ std::vector IC_SSHICM(const std::vector& d, // } // } // -// double p_value = static_cast(greater_count) / num_permutations; +// double p_value = static_cast(greater_count) / permutation_number; // // // Return a vector containing the true IC and p-value // return {true_IC, p_value}; diff --git a/src/IN_SSH.cpp b/src/IN_SSH.cpp index 50d19d2..dd4d34b 100644 --- a/src/IN_SSH.cpp +++ b/src/IN_SSH.cpp @@ -135,23 +135,26 @@ std::vector IN_SSHICM(const std::vector& d, // Step 2: Generate random permutations of d and compute IN_SSH for each std::vector IN_SSH_results(permutation_number, 0.0); // Store IN_SSH values for each permutation - // Initialize the random number generator with the seed - std::mt19937 gen(seed); + // Step 3: Generate a random seed using the input seed + std::mt19937 seed_gen(seed); // Initialize random number generator with the input seed + // std::uniform_int_distribution<> dis(1, std::numeric_limits::max()); // Define a distribution for random integers + std::uniform_int_distribution<> dis(1, 100); + int randomseed = dis(seed_gen); // Generate a random integer using the input seed - // Step 3: Perform parallel computation + // Step 4: Perform parallel computation RcppThread::parallelFor(0, permutation_number, [&](size_t i) { - // Generate a unique seed for each permutation - std::mt19937 local_gen(seed + i); // Modify seed for each thread + // Step 4.1: Generate a unique seed for each permutation by adding the iteration index to the randomseed + std::mt19937 local_gen(randomseed + i); // Modify seed for each thread - // Step 3.1: Permute d + // Step 4.2: Permute d std::vector permuted_d = d; // Copy the original d std::shuffle(permuted_d.begin(), permuted_d.end(), local_gen); // Shuffle based on the unique seed for each thread - // Step 3.2: Compute IN_SSH for the permuted d + // Step 4.3: Compute IN_SSH for the permuted d IN_SSH_results[i] = IN_SSH(permuted_d, s); }); - // Step 4: Compute p-value by comparing permuted IN_SSH values to the true IN_SSH value + // Step 5: Compute p-value by comparing permuted IN_SSH values to the true IN_SSH value int greater_count = 0; for (size_t i = 0; i < IN_SSH_results.size(); ++i) { if (IN_SSH_results[i] >= true_IN_SSH) { @@ -164,3 +167,49 @@ std::vector IN_SSHICM(const std::vector& d, // Return a vector containing the true IN_SSH and p-value return {true_IN_SSH, p_value}; } + +// // IN_SSHICM: Parallel computation of IN_SSH over permutations, returning IN_SSH value and p-value +// // [[Rcpp::export]] +// std::vector IN_SSHICM(const std::vector& d, +// const std::vector& s, +// unsigned int seed, +// int permutation_number) { +// if (s.size() != d.size()) { +// throw std::invalid_argument("Vectors s and d must have the same length."); +// } +// +// // Step 1: Calculate the true IN_SSH value using the original d and s +// double true_IN_SSH = IN_SSH(d, s); +// +// // Step 2: Generate random permutations of d and compute IN_SSH for each +// std::vector IN_SSH_results(permutation_number, 0.0); // Store IN_SSH values for each permutation +// +// // Initialize the random number generator with the seed +// std::mt19937 gen(seed); +// +// // Step 3: Perform parallel computation +// RcppThread::parallelFor(0, permutation_number, [&](size_t i) { +// // Generate a unique seed for each permutation +// std::mt19937 local_gen(seed + i); // Modify seed for each thread +// +// // Step 3.1: Permute d +// std::vector permuted_d = d; // Copy the original d +// std::shuffle(permuted_d.begin(), permuted_d.end(), local_gen); // Shuffle based on the unique seed for each thread +// +// // Step 3.2: Compute IN_SSH for the permuted d +// IN_SSH_results[i] = IN_SSH(permuted_d, s); +// }); +// +// // Step 4: Compute p-value by comparing permuted IN_SSH values to the true IN_SSH value +// int greater_count = 0; +// for (size_t i = 0; i < IN_SSH_results.size(); ++i) { +// if (IN_SSH_results[i] >= true_IN_SSH) { +// greater_count++; +// } +// } +// +// double p_value = static_cast(greater_count) / permutation_number; +// +// // Return a vector containing the true IN_SSH and p-value +// return {true_IN_SSH, p_value}; +// }