From eeda105ef06dc25d6552face779747f5b3764aae Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 01/26] enhance ccd methods --- inst/include/coupling/ccd.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/inst/include/coupling/ccd.hpp b/inst/include/coupling/ccd.hpp index e605273..20401c3 100644 --- a/inst/include/coupling/ccd.hpp +++ b/inst/include/coupling/ccd.hpp @@ -233,7 +233,7 @@ inline double ccd_c_single( double denom = (p - 1) * p / 2.0; double term1 = 1.0 - (sum_dist / denom); - if (term1 < 0) term1 = 0; + // if (term1 < 0) term1 = 0; double max_u = *std::max_element(vec.begin(), vec.end()); From a1979bbbb2396f67cda49abcb2b16e2a589c0cc3 Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 02/26] draft metacoupling method --- inst/include/coupling/ccd.hpp | 1 + inst/include/coupling/metacoupling.hpp | 170 +++++++++++++++++++++++++ 2 files changed, 171 insertions(+) create mode 100644 inst/include/coupling/metacoupling.hpp diff --git a/inst/include/coupling/ccd.hpp b/inst/include/coupling/ccd.hpp index 20401c3..c3ff563 100644 --- a/inst/include/coupling/ccd.hpp +++ b/inst/include/coupling/ccd.hpp @@ -178,6 +178,7 @@ #include #include #include +#include "coupling/numericutils.hpp" namespace coupling { diff --git a/inst/include/coupling/metacoupling.hpp b/inst/include/coupling/metacoupling.hpp new file mode 100644 index 0000000..125f85a --- /dev/null +++ b/inst/include/coupling/metacoupling.hpp @@ -0,0 +1,170 @@ +#ifndef COUPLING_METACOUPLING_HPP +#define COUPLING_METACOUPLING_HPP + +#include +#include +#include +#include +#include +#include +#include +#include "coupling/ccd.hpp" + +namespace coupling +{ + +namespace metacoupling +{ + +// mean +inline double mean(const std::vector& v) { + double s = std::accumulate(v.begin(), v.end(), 0.0); + return s / v.size(); +} + +inline double ccd_c_single( + const std::vector& vec, + const std::string& method = "standard" +) { + size_t p = vec.size(); // number of U values + + double C_val = 0.0; + + // ========================= + // standard + // ========================= + if (method == "standard") { + double prod_sum = 1.0; + for (double u : vec) { + // if (u <= 0) throw std::runtime_error("Values must be positive."); + prod_sum *= u; + } + + double geo_mean = std::pow(prod_sum, 1.0/p); + double arith_mean = mean(vec); + + if (coupling::numericutils::doubleNearlyEqual(arith_mean, 0.0)) { + return 0.0; + } + + C_val = geo_mean / arith_mean; + } + + // ========================= + // wang + // ========================= + else if (method == "wang") { + double sum_dist = 0.0; + + for (size_t i = 0; i < p - 1; ++i) { + for (size_t j = i + 1; j < p; ++j) { + sum_dist += std::abs(vec[i] - vec[j]); + } + } + + double denom = (p - 1) * p / 2.0; + double term1 = 1.0 - (sum_dist / denom); + // if (term1 < 0) term1 = 0; + + double max_u = *std::max_element(vec.begin(), vec.end()); + + if (coupling::numericutils::doubleNearlyEqual(max_u, 0.0)) { + return 0.0; + } + + double prod = 1.0; + for (double u : vec) { + prod *= (u / max_u); + } + + double term2 = std::pow(prod, 1.0 / (p - 1)); + + C_val = std::sqrt(term1 * term2); + } + + // ========================= + // fan + // ========================= + else if (method == "fan") { + double sum_u = std::accumulate(vec.begin(), vec.end(), 0.0); + // if (coupling::numericutils::doubleNearlyEqual(sum_u, 0.0)) { + // return 0.0; + // } + + double sum_u2 = 0.0; + for (double u : vec) { + sum_u2 += u * u; + } + // if (coupling::numericutils::doubleNearlyEqual(sum_u2, 0.0)) { + // return 0.0; + // } + + double numerator = p * sum_u2 - sum_u * sum_u; + double denom = p * p; + + double val = numerator / denom; + if (val < 0) val = 0; + + C_val = 1.0 - 2.0 * std::sqrt(val); + } + + else { + throw std::invalid_argument("Unknown method"); + } + + return std::clamp(C_val, 0.0, 1.0); +} + +inline std::vector ccd_c( + const std::vector>& mat, + const std::string& method = "standard", + size_t threads = 1 +) { + size_t n_units = mat.size(); + if (n_units == 0) return {}; + + std::vector result(n_units, 0.0); + + if (threads <= 1) { + for (size_t i = 0; i < n_units; ++i) { + result[i] = ccd_c_single(mat[i], method); + } + } else { + RcppThread::parallelFor(0, n_units, [&](size_t i) { + result[i] = ccd_c_single(mat[i], method); + }, threads); + } + + return result; +} + +inline std::vector> ccd( + const std::vector>& mat, + const std::vector& weight, + const std::string& method = "standard", + size_t threads = 1 +) { + size_t n_units = mat.size(); // number of unit + size_t p = mat[0].size(); // number of U values per unit + + std::vector> result(2, std::vector(n_units, 0.0)); + + std::vector C_vals = ccd_c(mat, method, threads); // C values + result[0] = C_vals; + + for (size_t i = 0; i < n_units; ++i) { + double T_val = 0.0; + for (size_t j = 0; j < p; ++j) { + T_val += weight[j] * mat[i][j]; // T = \sum_{j=1}^n w_j \times U_j + } + result[1][i] = std::sqrt(C_vals[i] * T_val); + } + + return result; +} + +} // namespace ccd + +} + +#endif // COUPLING_METACOUPLING_HPP From 72dc84f79917553f89610538489f703ffaa5287c Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 03/26] draft metacoupling method --- inst/include/coupling/metacoupling.hpp | 108 ++----------------------- 1 file changed, 6 insertions(+), 102 deletions(-) diff --git a/inst/include/coupling/metacoupling.hpp b/inst/include/coupling/metacoupling.hpp index 125f85a..6e17a3d 100644 --- a/inst/include/coupling/metacoupling.hpp +++ b/inst/include/coupling/metacoupling.hpp @@ -4,6 +4,7 @@ #include #include #include +#include #include #include #include @@ -16,114 +17,17 @@ namespace coupling namespace metacoupling { -// mean -inline double mean(const std::vector& v) { - double s = std::accumulate(v.begin(), v.end(), 0.0); - return s / v.size(); -} - -inline double ccd_c_single( - const std::vector& vec, - const std::string& method = "standard" -) { - size_t p = vec.size(); // number of U values - - double C_val = 0.0; - - // ========================= - // standard - // ========================= - if (method == "standard") { - double prod_sum = 1.0; - for (double u : vec) { - // if (u <= 0) throw std::runtime_error("Values must be positive."); - prod_sum *= u; - } - - double geo_mean = std::pow(prod_sum, 1.0/p); - double arith_mean = mean(vec); - - if (coupling::numericutils::doubleNearlyEqual(arith_mean, 0.0)) { - return 0.0; - } - - C_val = geo_mean / arith_mean; - } - - // ========================= - // wang - // ========================= - else if (method == "wang") { - double sum_dist = 0.0; - - for (size_t i = 0; i < p - 1; ++i) { - for (size_t j = i + 1; j < p; ++j) { - sum_dist += std::abs(vec[i] - vec[j]); - } - } - - double denom = (p - 1) * p / 2.0; - double term1 = 1.0 - (sum_dist / denom); - // if (term1 < 0) term1 = 0; - - double max_u = *std::max_element(vec.begin(), vec.end()); - - if (coupling::numericutils::doubleNearlyEqual(max_u, 0.0)) { - return 0.0; - } - - double prod = 1.0; - for (double u : vec) { - prod *= (u / max_u); - } - - double term2 = std::pow(prod, 1.0 / (p - 1)); - - C_val = std::sqrt(term1 * term2); - } - - // ========================= - // fan - // ========================= - else if (method == "fan") { - double sum_u = std::accumulate(vec.begin(), vec.end(), 0.0); - // if (coupling::numericutils::doubleNearlyEqual(sum_u, 0.0)) { - // return 0.0; - // } - - double sum_u2 = 0.0; - for (double u : vec) { - sum_u2 += u * u; - } - // if (coupling::numericutils::doubleNearlyEqual(sum_u2, 0.0)) { - // return 0.0; - // } - - double numerator = p * sum_u2 - sum_u * sum_u; - double denom = p * p; - - double val = numerator / denom; - if (val < 0) val = 0; - - C_val = 1.0 - 2.0 * std::sqrt(val); - } - - else { - throw std::invalid_argument("Unknown method"); - } - - return std::clamp(C_val, 0.0, 1.0); -} - -inline std::vector ccd_c( +inline std::vector metacoupling_c( const std::vector>& mat, + const std::vector>& swm_peri, + const std::vector>& swm_tele, const std::string& method = "standard", size_t threads = 1 ) { size_t n_units = mat.size(); - if (n_units == 0) return {}; - std::vector result(n_units, 0.0); + std::vector result(3, std::numeric_limits::quiet_NaN()); + if (n_units == 0) return result; if (threads <= 1) { for (size_t i = 0; i < n_units; ++i) { From e5ebf119eb3016202a3665360c6326c450fd0c0b Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 04/26] draft metacoupling method --- inst/include/coupling/metacoupling.hpp | 42 ++++++-------------------- 1 file changed, 9 insertions(+), 33 deletions(-) diff --git a/inst/include/coupling/metacoupling.hpp b/inst/include/coupling/metacoupling.hpp index 6e17a3d..f063cde 100644 --- a/inst/include/coupling/metacoupling.hpp +++ b/inst/include/coupling/metacoupling.hpp @@ -17,7 +17,7 @@ namespace coupling namespace metacoupling { -inline std::vector metacoupling_c( +inline std::vector> metacoupling_c( const std::vector>& mat, const std::vector>& swm_peri, const std::vector>& swm_tele, @@ -26,48 +26,24 @@ inline std::vector metacoupling_c( ) { size_t n_units = mat.size(); - std::vector result(3, std::numeric_limits::quiet_NaN()); - if (n_units == 0) return result; + std::vector> result( + n_units, + std::vector(3, std::numeric_limits::quiet_NaN())); if (threads <= 1) { - for (size_t i = 0; i < n_units; ++i) { - result[i] = ccd_c_single(mat[i], method); + for (size_t s = 0; s < n_units; ++s) { + result[s][0] = coupling::ccd::ccd_c_single(mat[s], method); } } else { - RcppThread::parallelFor(0, n_units, [&](size_t i) { - result[i] = ccd_c_single(mat[i], method); + RcppThread::parallelFor(0, n_units, [&](size_t s) { + result[s][0] = coupling::ccd::ccd_c_single(mat[s], method); }, threads); } return result; } -inline std::vector> ccd( - const std::vector>& mat, - const std::vector& weight, - const std::string& method = "standard", - size_t threads = 1 -) { - size_t n_units = mat.size(); // number of unit - size_t p = mat[0].size(); // number of U values per unit - - std::vector> result(2, std::vector(n_units, 0.0)); - - std::vector C_vals = ccd_c(mat, method, threads); // C values - result[0] = C_vals; - - for (size_t i = 0; i < n_units; ++i) { - double T_val = 0.0; - for (size_t j = 0; j < p; ++j) { - T_val += weight[j] * mat[i][j]; // T = \sum_{j=1}^n w_j \times U_j - } - result[1][i] = std::sqrt(C_vals[i] * T_val); - } - - return result; -} - -} // namespace ccd +} // namespace metacoupling } From 8d8cf1abd78d7cf7369069089fb3c1ba9c14c433 Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 05/26] draft metacoupling method --- inst/include/coupling/metacoupling.hpp | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/inst/include/coupling/metacoupling.hpp b/inst/include/coupling/metacoupling.hpp index f063cde..df4a4ac 100644 --- a/inst/include/coupling/metacoupling.hpp +++ b/inst/include/coupling/metacoupling.hpp @@ -24,7 +24,8 @@ inline std::vector> metacoupling_c( const std::string& method = "standard", size_t threads = 1 ) { - size_t n_units = mat.size(); + size_t n_units = mat.size(); // number of observationa + size_t p = mat[0].size(); // number of U values std::vector> result( n_units, From 9c62f0f9d58c6128223691fad18bb0082d0d201c Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 06/26] rename rcpp wrappers file --- src/{CCD.cpp => wrappers.cpp} | 0 1 file changed, 0 insertions(+), 0 deletions(-) rename src/{CCD.cpp => wrappers.cpp} (100%) diff --git a/src/CCD.cpp b/src/wrappers.cpp similarity index 100% rename from src/CCD.cpp rename to src/wrappers.cpp From 356fc9368b35f20632a7cabf989cc712a5dd21b9 Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 07/26] draft metacoupling method --- inst/include/coupling/metacoupling.hpp | 130 +++++++++++++++++++++++-- 1 file changed, 121 insertions(+), 9 deletions(-) diff --git a/inst/include/coupling/metacoupling.hpp b/inst/include/coupling/metacoupling.hpp index df4a4ac..6effe1e 100644 --- a/inst/include/coupling/metacoupling.hpp +++ b/inst/include/coupling/metacoupling.hpp @@ -17,28 +17,140 @@ namespace coupling namespace metacoupling { -inline std::vector> metacoupling_c( +inline std::vector> metacoupling_c( const std::vector>& mat, const std::vector>& swm_peri, const std::vector>& swm_tele, const std::string& method = "standard", size_t threads = 1 ) { - size_t n_units = mat.size(); // number of observationa - size_t p = mat[0].size(); // number of U values + // Number of spatial units (observations) + size_t n_units = mat.size(); + // Number of indicators (system dimensions) + size_t p = mat[0].size(); + + // Output: + // [i][0] = intra-coupling (original CCD for unit i) + // [i][1] = peri-coupling + // [i][2] = tele-coupling std::vector> result( n_units, - std::vector(3, std::numeric_limits::quiet_NaN())); - + std::vector(3, std::numeric_limits::quiet_NaN()) + ); + + // Total number of full binary combinations = 2^p + // Each bit represents whether a dimension comes from unit i or unit j + size_t full_perm = 1ULL << p; + + // ============================================================ + // Compute coupling between a pair (i, j) + // using HALF of the full permutation space: 2^(p-1) + // ============================================================ + auto compute_pair = [&](size_t i, size_t j) { + + double sum_c = 0.0; + size_t count = 0; + + // Iterate over all 2^p binary masks + for (size_t mask = 0; mask < full_perm; ++mask) { + + // ---------------------------------------------------- + // Symmetry reduction: + // For any mask, its complement (~mask) produces + // the same combination but with i and j swapped. + // + // Example: + // mask = 0101 + // complement= 1010 + // + // These two correspond to equivalent mixed systems. + // + // To avoid double counting, we only keep one of them: + // enforce: mask <= complement + // ---------------------------------------------------- + size_t complement = (~mask) & (full_perm - 1); + + if (mask > complement) continue; + + // Construct a mixed indicator vector + std::vector vec(p); + + // ---------------------------------------------------- + // Bitwise construction: + // For each dimension k: + // + // if bit k == 1 → take value from unit i + // if bit k == 0 → take value from unit j + // + // This ensures: + // - all possible combinations are explored + // - no dimension is artificially fixed + // ---------------------------------------------------- + for (size_t k = 0; k < p; ++k) { + if (mask & (1ULL << k)) { + vec[k] = mat[i][k]; + } else { + vec[k] = mat[j][k]; + } + } + + // Compute CCD coupling for this mixed system + sum_c += coupling::ccd::ccd_c_single(vec, method); + ++count; + } + + // The number of retained combinations is exactly 2^(p-1) + return sum_c / static_cast(count); + }; + + // ============================================================ + // Main worker for each spatial unit s + // ============================================================ + auto worker = [&](size_t s) { + + // 1. Intra-coupling (original system) + result[s][0] = coupling::ccd::ccd_c_single(mat[s], method); + + double peri_sum = 0.0; + double tele_sum = 0.0; + + // 2. Iterate over all other units + for (size_t j = 0; j < n_units; ++j) { + + double w_peri = swm_peri[s][j]; + double w_tele = swm_tele[s][j]; + + // Skip if no spatial interaction + if (w_peri == 0.0 && w_tele == 0.0) continue; + + // Compute symmetric coupling between (s, j) + double c_val = compute_pair(s, j); + + // 3. Weighted aggregation + // Since SWM is symmetric, no need to compute (j, s) separately + if (w_peri != 0.0) { + peri_sum += w_peri * c_val; + } + + if (w_tele != 0.0) { + tele_sum += w_tele * c_val; + } + } + + result[s][1] = peri_sum; + result[s][2] = tele_sum; + }; + + // ============================================================ + // Parallel or serial execution + // ============================================================ if (threads <= 1) { for (size_t s = 0; s < n_units; ++s) { - result[s][0] = coupling::ccd::ccd_c_single(mat[s], method); + worker(s); } } else { - RcppThread::parallelFor(0, n_units, [&](size_t s) { - result[s][0] = coupling::ccd::ccd_c_single(mat[s], method); - }, threads); + RcppThread::parallelFor(0, n_units, worker, threads); } return result; From e8abd49769898462b3099f4120a9fa7ceb2b23c2 Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 08/26] draft metacoupling method --- inst/include/coupling/metacoupling.hpp | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/inst/include/coupling/metacoupling.hpp b/inst/include/coupling/metacoupling.hpp index 6effe1e..45e59b9 100644 --- a/inst/include/coupling/metacoupling.hpp +++ b/inst/include/coupling/metacoupling.hpp @@ -121,19 +121,20 @@ inline std::vector> metacoupling_c( double w_peri = swm_peri[s][j]; double w_tele = swm_tele[s][j]; - // Skip if no spatial interaction - if (w_peri == 0.0 && w_tele == 0.0) continue; + // Skip if no interaction + if (coupling::numericutils::doubleNearlyEqual(w_peri, 0.0) && + coupling::numericutils::doubleNearlyEqual(w_tele, 0.0)) continue; // Compute symmetric coupling between (s, j) double c_val = compute_pair(s, j); // 3. Weighted aggregation // Since SWM is symmetric, no need to compute (j, s) separately - if (w_peri != 0.0) { + if (!coupling::numericutils::doubleNearlyEqual(w_peri, 0.0)) { peri_sum += w_peri * c_val; } - if (w_tele != 0.0) { + if (!coupling::numericutils::doubleNearlyEqual(w_tele, 0.0)) { tele_sum += w_tele * c_val; } } From a0a3b0dc3f04d8a20bce3a9dd860bb16317afcd0 Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 09/26] draft metacoupling method --- inst/include/coupling/metacoupling.hpp | 55 +++++++++++--------------- 1 file changed, 23 insertions(+), 32 deletions(-) diff --git a/inst/include/coupling/metacoupling.hpp b/inst/include/coupling/metacoupling.hpp index 45e59b9..2f06435 100644 --- a/inst/include/coupling/metacoupling.hpp +++ b/inst/include/coupling/metacoupling.hpp @@ -39,54 +39,45 @@ inline std::vector> metacoupling_c( std::vector(3, std::numeric_limits::quiet_NaN()) ); - // Total number of full binary combinations = 2^p - // Each bit represents whether a dimension comes from unit i or unit j + // Total number of binary combinations = 2^p size_t full_perm = 1ULL << p; // ============================================================ - // Compute coupling between a pair (i, j) - // using HALF of the full permutation space: 2^(p-1) + // Compute metacoupling between unit i and j + // + // Definition: + // - Enumerate ALL mixed systems + // - Each dimension independently comes from i or j + // - EXCLUDE: + // mask == 0 → all from j + // mask == full_perm-1 → all from i + // + // Total combinations used = 2^p - 2 // ============================================================ auto compute_pair = [&](size_t i, size_t j) { double sum_c = 0.0; size_t count = 0; - // Iterate over all 2^p binary masks for (size_t mask = 0; mask < full_perm; ++mask) { // ---------------------------------------------------- - // Symmetry reduction: - // For any mask, its complement (~mask) produces - // the same combination but with i and j swapped. - // - // Example: - // mask = 0101 - // complement= 1010 - // - // These two correspond to equivalent mixed systems. - // - // To avoid double counting, we only keep one of them: - // enforce: mask <= complement + // Remove trivial (non-mixed) systems // ---------------------------------------------------- - size_t complement = (~mask) & (full_perm - 1); - - if (mask > complement) continue; - - // Construct a mixed indicator vector - std::vector vec(p); + if (mask == 0 || mask == (full_perm - 1)) continue; // ---------------------------------------------------- - // Bitwise construction: - // For each dimension k: + // Construct mixed system vector // - // if bit k == 1 → take value from unit i - // if bit k == 0 → take value from unit j + // Bit k of mask determines source: + // 1 → take from unit i + // 0 → take from unit j // - // This ensures: - // - all possible combinations are explored - // - no dimension is artificially fixed + // Example (p = 3): + // mask = 101 → [i, j, i] // ---------------------------------------------------- + std::vector vec(p); + for (size_t k = 0; k < p; ++k) { if (mask & (1ULL << k)) { vec[k] = mat[i][k]; @@ -95,12 +86,12 @@ inline std::vector> metacoupling_c( } } - // Compute CCD coupling for this mixed system + // Compute CCD coupling sum_c += coupling::ccd::ccd_c_single(vec, method); ++count; } - // The number of retained combinations is exactly 2^(p-1) + // count should be exactly 2^p - 2 return sum_c / static_cast(count); }; From 7ca0a3740ccb664dd89ad8ca850903503ab407a9 Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 10/26] draft metacoupling method --- inst/include/coupling/ccd.hpp | 1 + inst/include/coupling/metacoupling.hpp | 2 ++ 2 files changed, 3 insertions(+) diff --git a/inst/include/coupling/ccd.hpp b/inst/include/coupling/ccd.hpp index c3ff563..d978cb8 100644 --- a/inst/include/coupling/ccd.hpp +++ b/inst/include/coupling/ccd.hpp @@ -197,6 +197,7 @@ inline double ccd_c_single( const std::string& method = "standard" ) { size_t p = vec.size(); // number of U values + if (p <= 1) return std::numeric_limits::quiet_NaN(); double C_val = 0.0; diff --git a/inst/include/coupling/metacoupling.hpp b/inst/include/coupling/metacoupling.hpp index 2f06435..110e5e5 100644 --- a/inst/include/coupling/metacoupling.hpp +++ b/inst/include/coupling/metacoupling.hpp @@ -39,6 +39,8 @@ inline std::vector> metacoupling_c( std::vector(3, std::numeric_limits::quiet_NaN()) ); + if (p <= 1) return result; + // Total number of binary combinations = 2^p size_t full_perm = 1ULL << p; From 8acdae15e1e342c4140390c2ab6cd35603a0836b Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 11/26] draft metacoupling method --- inst/include/coupling/metacoupling.hpp | 110 +++++++++++++++++++++++++ 1 file changed, 110 insertions(+) diff --git a/inst/include/coupling/metacoupling.hpp b/inst/include/coupling/metacoupling.hpp index 110e5e5..d989539 100644 --- a/inst/include/coupling/metacoupling.hpp +++ b/inst/include/coupling/metacoupling.hpp @@ -150,6 +150,116 @@ inline std::vector> metacoupling_c( return result; } +inline std::vector> metacoupling( + const std::vector>& mat, + const std::vector>& swm_peri, + const std::vector>& swm_tele, + const std::vector& weight, + const std::string& method = "standard", + size_t threads = 1 +) { + size_t n_units = mat.size(); + size_t p = mat[0].size(); + + std::vector> result( + n_units, + std::vector(6, std::numeric_limits::quiet_NaN()) + ); + + if (p <= 1) return result; + + size_t full_perm = 1ULL << p; + + // ============================================================ + // Step 1: reuse existing C computation + // ============================================================ + auto C_res = metacoupling_c(mat, swm_peri, swm_tele, method, threads); + + // ============================================================ + // Step 2: T computation + // ============================================================ + auto compute_T_intra = [&](size_t i) { + double t = 0.0; + for (size_t k = 0; k < p; ++k) { + t += weight[k] * mat[i][k]; + } + return t; + }; + + auto compute_pair_T = [&](size_t i, size_t j) { + + double sum_t = 0.0; + size_t count = 0; + + for (size_t mask = 0; mask < full_perm; ++mask) { + + if (mask == 0 || mask == (full_perm - 1)) continue; + + double t_val = 0.0; + + for (size_t k = 0; k < p; ++k) { + if (mask & (1ULL << k)) { + t_val += weight[k] * mat[i][k]; + } else { + t_val += weight[k] * mat[j][k]; + } + } + + sum_t += t_val; + ++count; + } + + return sum_t / static_cast(count); + }; + + // ============================================================ + // Step 3: assemble results + // ============================================================ + for (size_t s = 0; s < n_units; ++s) { + + // ---- intra ---- + double C_intra = C_res[s][0]; + double T_intra = compute_T_intra(s); + + result[s][0] = C_intra; + result[s][1] = std::sqrt(C_intra * T_intra); + + double peri_T_sum = 0.0; + double tele_T_sum = 0.0; + + // ---- interactions ---- + for (size_t j = 0; j < n_units; ++j) { + + double w_peri = swm_peri[s][j]; + double w_tele = swm_tele[s][j]; + + if (coupling::numericutils::doubleNearlyEqual(w_peri, 0.0) && + coupling::numericutils::doubleNearlyEqual(w_tele, 0.0)) continue; + + double t_val = compute_pair_T(s, j); + + if (!coupling::numericutils::doubleNearlyEqual(w_peri, 0.0)) { + peri_T_sum += w_peri * t_val; + } + + if (!coupling::numericutils::doubleNearlyEqual(w_tele, 0.0)) { + tele_T_sum += w_tele * t_val; + } + } + + double C_peri = C_res[s][1]; + double C_tele = C_res[s][2]; + + result[s][2] = C_peri; + result[s][3] = std::sqrt(C_peri * peri_T_sum); + + result[s][4] = C_tele; + result[s][5] = std::sqrt(C_tele * tele_T_sum); + } + + return result; +} + } // namespace metacoupling } From b1004feb0f4b322ae9992c791157e76bac764bae Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 12/26] draft metacoupling method --- inst/include/coupling/metacoupling.hpp | 1 + 1 file changed, 1 insertion(+) diff --git a/inst/include/coupling/metacoupling.hpp b/inst/include/coupling/metacoupling.hpp index d989539..1c0dc1c 100644 --- a/inst/include/coupling/metacoupling.hpp +++ b/inst/include/coupling/metacoupling.hpp @@ -10,6 +10,7 @@ #include #include #include "coupling/ccd.hpp" +#include "coupling/numericutils.hpp" namespace coupling { From 3f8884f4b453d7ad337f9db43d181817fc4f252b Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 13/26] draft metacoupling method --- inst/include/coupling/metacoupling.hpp | 167 +++++++++++++++++++++++++ 1 file changed, 167 insertions(+) diff --git a/inst/include/coupling/metacoupling.hpp b/inst/include/coupling/metacoupling.hpp index 1c0dc1c..7e00f20 100644 --- a/inst/include/coupling/metacoupling.hpp +++ b/inst/include/coupling/metacoupling.hpp @@ -1,3 +1,170 @@ +/****************************************************************************** + * File: metacoupling.hpp + * + * Metacoupling Model (Extension of CCD) + * ------------------------------------------------ + * + * This module implements the metacoupling framework, an extension of the + * Coupling Coordination Degree (CCD) model, designed to measure both + * intra-system and inter-system coupling relationships across spatial units. + * + * The metacoupling model integrates: + * + * • intra-unit coupling (within a spatial unit) + * • peri-coupling (short-distance interactions) + * • tele-coupling (long-distance interactions) + * + * It further extends the CCD framework by introducing a permutation-based + * cross-system coupling mechanism. + * + * --------------------------------------------------------------------------- + * Mathematical formulation + * --------------------------------------------------------------------------- + * + * Let there be: + * + * • n spatial units + * • p indicators per unit + * + * For unit i: + * + * U_i = {U_{i1}, U_{i2}, ..., U_{ip}} + * + * --------------------------------------------------------------------------- + * 1. Intra-system coupling (C_intra) + * --------------------------------------------------------------------------- + * + * Standard CCD coupling computed for each unit: + * + * C_intra(i) = CCD(U_i) + * + * where CCD is defined in ccd.hpp. + * + * --------------------------------------------------------------------------- + * 2. Inter-system metacoupling (C_ij) + * --------------------------------------------------------------------------- + * + * For two units i and j, metacoupling is defined by constructing + * mixed systems through all possible indicator-level combinations: + * + * V(mask) = {v_1, v_2, ..., v_p} + * + * where: + * + * v_k = U_{ik} if bit k in mask = 1 + * v_k = U_{jk} if bit k in mask = 0 + * + * All binary masks are enumerated: + * + * mask ∈ {0, 1}^p + * + * EXCLUDING: + * + * • mask = 0 (all indicators from j) + * • mask = 2^p − 1 (all indicators from i) + * + * Total combinations used: + * + * 2^p − 2 + * + * The inter-system coupling is: + * + * C_ij = mean over all valid masks of CCD(V(mask)) + * + * --------------------------------------------------------------------------- + * 3. Spatial aggregation (peri / tele) + * --------------------------------------------------------------------------- + * + * Given spatial weight matrices: + * + * W_peri (short-distance interactions) + * W_tele (long-distance interactions) + * + * Both are assumed to be symmetric. + * + * Aggregated coupling: + * + * C_peri(i) = Σ_j W_peri(i,j) * C_ij + * C_tele(i) = Σ_j W_tele(i,j) * C_ij + * + * --------------------------------------------------------------------------- + * 4. Composite development index (T) + * --------------------------------------------------------------------------- + * + * For each system: + * + * T = Σ_k w_k * U_k + * + * For inter-system: + * + * T_ij = mean over all permutations of T(V(mask)) + * + * --------------------------------------------------------------------------- + * 5. Coupling coordination degree (D) + * --------------------------------------------------------------------------- + * + * Final metric: + * + * D = sqrt(C × T) + * + * Applied to: + * + * • intra-system + * • peri-coupling + * • tele-coupling + * + * --------------------------------------------------------------------------- + * Input data format + * --------------------------------------------------------------------------- + * + * mat : vector> + * + * • rows = spatial units + * • columns = indicators + * • values MUST be normalized to [0,1] + * + * swm_peri / swm_tele: + * + * • symmetric spatial weight matrices (n × n) + * + * weight: + * + * • length = number of indicators + * + * --------------------------------------------------------------------------- + * Output + * --------------------------------------------------------------------------- + * + * metacoupling_c: + * + * result[i] = {C_intra, C_peri, C_tele} + * + * metacoupling: + * + * result[i] = { + * C_intra, D_intra, + * C_peri, D_peri, + * C_tele, D_tele + * } + * + * --------------------------------------------------------------------------- + * Notes + * --------------------------------------------------------------------------- + * + * • The number of inter-system combinations grows exponentially (2^p − 2) + * + * • This framework captures cross-unit structural interactions + * at the indicator level, rather than treating systems as indivisible units + * + * • Input normalization is critical for meaningful interpretation + * + * • Computational cost increases rapidly with p + * + * --------------------------------------------------------------------------- + * Author: Wenbo Lyu (Github: @SpatLyu) + * License: GPL-3 + ******************************************************************************/ + #ifndef COUPLING_METACOUPLING_HPP #define COUPLING_METACOUPLING_HPP From 6d599bcbed83a2d9977ae245bfa70eb15a96761a Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 14/26] draft metacoupling method --- inst/include/coupling.h | 1 + 1 file changed, 1 insertion(+) diff --git a/inst/include/coupling.h b/inst/include/coupling.h index 0456ec5..6d73838 100644 --- a/inst/include/coupling.h +++ b/inst/include/coupling.h @@ -39,6 +39,7 @@ #include "coupling/numericutils.hpp" #include "coupling/ccd.hpp" +#include "coupling/metacoupling.hpp" // ============================================================ // Convenience Converters (Inline helpers for R/C++ interop) From 18f3a7c4db26080e9ab7b43fe8fe3b7c2648b82f Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 15/26] implement metacoupling wrappers --- src/wrappers.cpp | 30 ++++++++++++++++++++++++++++++ 1 file changed, 30 insertions(+) diff --git a/src/wrappers.cpp b/src/wrappers.cpp index 309c353..9ba950a 100644 --- a/src/wrappers.cpp +++ b/src/wrappers.cpp @@ -25,3 +25,33 @@ Rcpp::DataFrame RcppCCD(const Rcpp::NumericMatrix& mat, Rcpp::Named("D") = res[1] ); } + +// Wrapper function to preform meta-coupling analysis +// [[Rcpp::export(rng = false)]] +Rcpp::DataFrame RcppMetaCoupling( + const Rcpp::NumericMatrix& mat, + const Rcpp::NumericMatrix& swm_peri, + const Rcpp::NumericMatrix& swm_tele, + const Rcpp::NumericVector& weight, + const std::string& method = "standard", + int threads = 1) +{ + std::vector> mat_std = coupling::convert::mat_r2std(mat, true); + std::vector> swm_peri_std = coupling::convert::mat_r2std(swm_peri, true); + std::vector> swm_tele_std = coupling::convert::mat_r2std(swm_tele, true); + std::vector wt_std = Rcpp::as>(weight); + + // Call metacoupling function + std::vector> res = coupling::metacoupling::metacoupling( + mat_std, swm_peri_std, swm_tele_std, wt_std, method, static_cast(std::abs(threads)) + ); + + return Rcpp::DataFrame::create( + Rcpp::Named("Intra_C") = res[0], + Rcpp::Named("Peri_C") = res[1], + Rcpp::Named("Tele_C") = res[2], + Rcpp::Named("Intra_D") = res[3], + Rcpp::Named("Peri_D") = res[4], + Rcpp::Named("Tele_D") = res[5] + ); +} From 106f781df24c28947a9c834adaa34bce143f23a6 Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 16/26] implement metacoupling wrappers --- src/wrappers.cpp | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/wrappers.cpp b/src/wrappers.cpp index 9ba950a..6173dcd 100644 --- a/src/wrappers.cpp +++ b/src/wrappers.cpp @@ -48,10 +48,10 @@ Rcpp::DataFrame RcppMetaCoupling( return Rcpp::DataFrame::create( Rcpp::Named("Intra_C") = res[0], - Rcpp::Named("Peri_C") = res[1], - Rcpp::Named("Tele_C") = res[2], - Rcpp::Named("Intra_D") = res[3], - Rcpp::Named("Peri_D") = res[4], + Rcpp::Named("Intra_D") = res[1], + Rcpp::Named("Peri_C") = res[2], + Rcpp::Named("Peri_D") = res[3], + Rcpp::Named("Tele_C") = res[4], Rcpp::Named("Tele_D") = res[5] ); } From 0abfd72987461bfa1e4af7859b6977c4680ba0ea Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 17/26] implement metacoupling wrappers --- inst/include/coupling/metacoupling.hpp | 50 +++++++++++++++----------- 1 file changed, 30 insertions(+), 20 deletions(-) diff --git a/inst/include/coupling/metacoupling.hpp b/inst/include/coupling/metacoupling.hpp index 7e00f20..0a5addf 100644 --- a/inst/include/coupling/metacoupling.hpp +++ b/inst/include/coupling/metacoupling.hpp @@ -199,12 +199,12 @@ inline std::vector> metacoupling_c( size_t p = mat[0].size(); // Output: - // [i][0] = intra-coupling (original CCD for unit i) - // [i][1] = peri-coupling - // [i][2] = tele-coupling + // [0][i] = intra-coupling (original CCD for unit i) + // [1][i] = peri-coupling + // [2][i] = tele-coupling std::vector> result( - n_units, - std::vector(3, std::numeric_limits::quiet_NaN()) + 3, + std::vector(n_units, std::numeric_limits::quiet_NaN()) ); if (p <= 1) return result; @@ -271,7 +271,7 @@ inline std::vector> metacoupling_c( auto worker = [&](size_t s) { // 1. Intra-coupling (original system) - result[s][0] = coupling::ccd::ccd_c_single(mat[s], method); + result[0][s] = coupling::ccd::ccd_c_single(mat[s], method); double peri_sum = 0.0; double tele_sum = 0.0; @@ -300,8 +300,8 @@ inline std::vector> metacoupling_c( } } - result[s][1] = peri_sum; - result[s][2] = tele_sum; + result[1][s] = peri_sum; + result[2][s] = tele_sum; }; // ============================================================ @@ -328,10 +328,20 @@ inline std::vector> metacoupling( ) { size_t n_units = mat.size(); size_t p = mat[0].size(); - + + // ============================================================ + // Output: 6 × n_units + // + // [0] = Intra_C + // [1] = Intra_D + // [2] = Peri_C + // [3] = Peri_D + // [4] = Tele_C + // [5] = Tele_D + // ============================================================ std::vector> result( - n_units, - std::vector(6, std::numeric_limits::quiet_NaN()) + 6, + std::vector(n_units, std::numeric_limits::quiet_NaN()) ); if (p <= 1) return result; @@ -386,11 +396,11 @@ inline std::vector> metacoupling( for (size_t s = 0; s < n_units; ++s) { // ---- intra ---- - double C_intra = C_res[s][0]; + double C_intra = C_res[0][s]; double T_intra = compute_T_intra(s); - result[s][0] = C_intra; - result[s][1] = std::sqrt(C_intra * T_intra); + result[0][s] = C_intra; + result[1][s] = std::sqrt(C_intra * T_intra); double peri_T_sum = 0.0; double tele_T_sum = 0.0; @@ -415,14 +425,14 @@ inline std::vector> metacoupling( } } - double C_peri = C_res[s][1]; - double C_tele = C_res[s][2]; + double C_peri = C_res[1][s]; + double C_tele = C_res[2][s]; - result[s][2] = C_peri; - result[s][3] = std::sqrt(C_peri * peri_T_sum); + result[2][s] = C_peri; + result[3][s] = std::sqrt(C_peri * peri_T_sum); - result[s][4] = C_tele; - result[s][5] = std::sqrt(C_tele * tele_T_sum); + result[4][s] = C_tele; + result[5][s] = std::sqrt(C_tele * tele_T_sum); } return result; From be0059e5804ac97d74636ba9a2b5223c3754e8d8 Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 18/26] implement metacoupling wrappers --- inst/include/coupling/metacoupling.hpp | 21 ++++++++++++++------- 1 file changed, 14 insertions(+), 7 deletions(-) diff --git a/inst/include/coupling/metacoupling.hpp b/inst/include/coupling/metacoupling.hpp index 0a5addf..5fdf3fd 100644 --- a/inst/include/coupling/metacoupling.hpp +++ b/inst/include/coupling/metacoupling.hpp @@ -131,21 +131,28 @@ * * • length = number of indicators * - * --------------------------------------------------------------------------- + * --------------------------------------------------------------------------- * Output * --------------------------------------------------------------------------- * * metacoupling_c: * - * result[i] = {C_intra, C_peri, C_tele} + * result is a 3 × n_units matrix: + * + * result[0][i] = C_intra(i) + * result[1][i] = C_peri(i) + * result[2][i] = C_tele(i) * * metacoupling: * - * result[i] = { - * C_intra, D_intra, - * C_peri, D_peri, - * C_tele, D_tele - * } + * result is a 6 × n_units matrix: + * + * result[0][i] = C_intra(i) + * result[1][i] = D_intra(i) + * result[2][i] = C_peri(i) + * result[3][i] = D_peri(i) + * result[4][i] = C_tele(i) + * result[5][i] = D_tele(i) * * --------------------------------------------------------------------------- * Notes From b28dca0a2b6b018332797496668400c3705c8f8b Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 19/26] implement metacoupling wrappers --- R/RcppExports.R | 4 ++++ src/RcppExports.cpp | 16 ++++++++++++++++ 2 files changed, 20 insertions(+) diff --git a/R/RcppExports.R b/R/RcppExports.R index bcf759d..5f41588 100644 --- a/R/RcppExports.R +++ b/R/RcppExports.R @@ -5,3 +5,7 @@ RcppCCD <- function(mat, weight, method = "standard", threads = 1L) { .Call(`_coupling_RcppCCD`, mat, weight, method, threads) } +RcppMetaCoupling <- function(mat, swm_peri, swm_tele, weight, method = "standard", threads = 1L) { + .Call(`_coupling_RcppMetaCoupling`, mat, swm_peri, swm_tele, weight, method, threads) +} + diff --git a/src/RcppExports.cpp b/src/RcppExports.cpp index afb6ae8..8dae5df 100644 --- a/src/RcppExports.cpp +++ b/src/RcppExports.cpp @@ -25,9 +25,25 @@ BEGIN_RCPP return rcpp_result_gen; END_RCPP } +// RcppMetaCoupling +Rcpp::DataFrame RcppMetaCoupling(const Rcpp::NumericMatrix& mat, const Rcpp::NumericMatrix& swm_peri, const Rcpp::NumericMatrix& swm_tele, const Rcpp::NumericVector& weight, const std::string& method, int threads); +RcppExport SEXP _coupling_RcppMetaCoupling(SEXP matSEXP, SEXP swm_periSEXP, SEXP swm_teleSEXP, SEXP weightSEXP, SEXP methodSEXP, SEXP threadsSEXP) { +BEGIN_RCPP + Rcpp::RObject rcpp_result_gen; + Rcpp::traits::input_parameter< const Rcpp::NumericMatrix& >::type mat(matSEXP); + Rcpp::traits::input_parameter< const Rcpp::NumericMatrix& >::type swm_peri(swm_periSEXP); + Rcpp::traits::input_parameter< const Rcpp::NumericMatrix& >::type swm_tele(swm_teleSEXP); + Rcpp::traits::input_parameter< const Rcpp::NumericVector& >::type weight(weightSEXP); + Rcpp::traits::input_parameter< const std::string& >::type method(methodSEXP); + Rcpp::traits::input_parameter< int >::type threads(threadsSEXP); + rcpp_result_gen = Rcpp::wrap(RcppMetaCoupling(mat, swm_peri, swm_tele, weight, method, threads)); + return rcpp_result_gen; +END_RCPP +} static const R_CallMethodDef CallEntries[] = { {"_coupling_RcppCCD", (DL_FUNC) &_coupling_RcppCCD, 4}, + {"_coupling_RcppMetaCoupling", (DL_FUNC) &_coupling_RcppMetaCoupling, 6}, {NULL, NULL, 0} }; From 0a929256d5f9d09eb1f8c7579559675e6ee0223f Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 20/26] implement metacoupling wrappers --- inst/include/coupling/metacoupling.hpp | 2 ++ 1 file changed, 2 insertions(+) diff --git a/inst/include/coupling/metacoupling.hpp b/inst/include/coupling/metacoupling.hpp index 5fdf3fd..5e0b639 100644 --- a/inst/include/coupling/metacoupling.hpp +++ b/inst/include/coupling/metacoupling.hpp @@ -21,6 +21,8 @@ * Mathematical formulation * --------------------------------------------------------------------------- * + * (See https://github.com/stscl/coupling/discussions/8 for a better reading experience.) + * * Let there be: * * • n spatial units From ca9e3c21811146406f6d627e905ea79c41717e6e Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Wed, 17 Jun 2026 00:00:01 +0800 Subject: [PATCH 21/26] implement metacoupling wrappers --- R/metacoupling.R | 9 +++++++++ 1 file changed, 9 insertions(+) create mode 100644 R/metacoupling.R diff --git a/R/metacoupling.R b/R/metacoupling.R new file mode 100644 index 0000000..791e47c --- /dev/null +++ b/R/metacoupling.R @@ -0,0 +1,9 @@ +metacoupling = \(data, swm_peri = NULL, swm_tele = NULL, weight = NULL, + method = c("standard", "wang", "fan"), threads = 1){ + mat = as.matrix(data) + method = match.arg(method) + if (is.null(weight)) weight = rep(1, times = ncol(mat)) / ncol(mat) + if (is.null(swm_peri)) swm_peri = matrix(0, nrow(data), nrow(data)) + if (is.null(swm_tele)) swm_tele = matrix(0, nrow(data), nrow(data)) + return(RcppMetaCoupling(mat, swm_peri, swm_tele, weight, method, threads)) +} From 0b629afcf04a2ecf14841a2eace137c8427e6981 Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Thu, 18 Jun 2026 00:00:01 +0800 Subject: [PATCH 22/26] implement metacoupling wrappers --- R/metacoupling.R | 35 +++++++++++++++++++++++++++++++++++ 1 file changed, 35 insertions(+) diff --git a/R/metacoupling.R b/R/metacoupling.R index 791e47c..05b7100 100644 --- a/R/metacoupling.R +++ b/R/metacoupling.R @@ -1,3 +1,38 @@ +#' Metacoupling Analysis +#' +#' @inheritParams ccd +#' @param swm_peri A numeric matrix representing the **peri (local) spatial weight matrix**. +#' Must be square with dimension equal to `nrow(data)`. If `NULL`, a zero matrix is used. +#' @param swm_tele A numeric matrix representing the **tele (long-distance) spatial weight matrix**. +#' Must be square with dimension equal to `nrow(data)`. If `NULL`, a zero matrix is used. +#' +#' @return A data.frame with: +#' \itemize{ +#' \item `Intra_C`: intra-system coupling degree +#' \item `Intra_D`: intra-system coordination degree +#' \item `Peri_C`: peri-coupling degree +#' \item `Peri_D`: peri coordination degree +#' \item `Tele_C`: tele-coupling degree +#' \item `Tele_D`: tele coordination degree +#' } +#' +#' @export +#' +#' @details +#' Full model definitions and formulas are available at: +#' \url{https://github.com/stscl/coupling/discussions/8} +#' +#' @note +#' Input values should be normalized to `[0, 1]`. Spatial weight matrices are +#' typically symmetric. +#' +#' @examples +#' set.seed(42) +#' mat = matrix(runif(20), nrow = 5) +#' swm1 = matrix(runif(25), 5, 5) +#' swm2 = matrix(runif(25), 5, 5) +#' coupling::metacoupling(mat, swm1, swm2) +#' metacoupling = \(data, swm_peri = NULL, swm_tele = NULL, weight = NULL, method = c("standard", "wang", "fan"), threads = 1){ mat = as.matrix(data) From 59aae45f4f2ba45210248928c8477b1d3c4c7dd3 Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Thu, 18 Jun 2026 00:00:01 +0800 Subject: [PATCH 23/26] implement metacoupling wrappers --- NAMESPACE | 1 + man/metacoupling.Rd | 62 +++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 63 insertions(+) create mode 100644 man/metacoupling.Rd diff --git a/NAMESPACE b/NAMESPACE index 7f5f347..bc23093 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -1,4 +1,5 @@ # Generated by roxygen2: do not edit by hand export(ccd) +export(metacoupling) useDynLib(coupling, .registration = TRUE) diff --git a/man/metacoupling.Rd b/man/metacoupling.Rd new file mode 100644 index 0000000..1330fb5 --- /dev/null +++ b/man/metacoupling.Rd @@ -0,0 +1,62 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/metacoupling.R +\name{metacoupling} +\alias{metacoupling} +\title{Metacoupling Analysis} +\usage{ +metacoupling( + data, + swm_peri = NULL, + swm_tele = NULL, + weight = NULL, + method = c("standard", "wang", "fan"), + threads = 1 +) +} +\arguments{ +\item{data}{A numeric matrix or data.frame. Rows are observations, +columns are indicators.} + +\item{swm_peri}{A numeric matrix representing the \strong{peri (local) spatial weight matrix}. +Must be square with dimension equal to \code{nrow(data)}. If \code{NULL}, a zero matrix is used.} + +\item{swm_tele}{A numeric matrix representing the \strong{tele (long-distance) spatial weight matrix}. +Must be square with dimension equal to \code{nrow(data)}. If \code{NULL}, a zero matrix is used.} + +\item{weight}{Numeric vector of indicator weights. Must have length equal +to \code{ncol(data)}. If \code{NULL}, equal weights are used.} + +\item{method}{Coupling model. One of \code{"standard"}, \code{"wang"}, or \code{"fan"}.} + +\item{threads}{Number of threads used in computation.} +} +\value{ +A data.frame with: +\itemize{ +\item \code{Intra_C}: intra-system coupling degree +\item \code{Intra_D}: intra-system coordination degree +\item \code{Peri_C}: peri-coupling degree +\item \code{Peri_D}: peri coordination degree +\item \code{Tele_C}: tele-coupling degree +\item \code{Tele_D}: tele coordination degree +} +} +\description{ +Metacoupling Analysis +} +\details{ +Full model definitions and formulas are available at: +\url{https://github.com/stscl/coupling/discussions/8} +} +\note{ +Input values should be normalized to \verb{[0, 1]}. Spatial weight matrices are +typically symmetric. +} +\examples{ +set.seed(42) +mat = matrix(runif(20), nrow = 5) +swm1 = matrix(runif(25), 5, 5) +swm2 = matrix(runif(25), 5, 5) +coupling::metacoupling(mat, swm1, swm2) + +} From 8cd0eb380e7281f7843505818c9d81fd2c4ffd8b Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Thu, 18 Jun 2026 00:00:01 +0800 Subject: [PATCH 24/26] implement metacoupling wrappers --- R/metacoupling.R | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/R/metacoupling.R b/R/metacoupling.R index 05b7100..9168d51 100644 --- a/R/metacoupling.R +++ b/R/metacoupling.R @@ -29,8 +29,8 @@ #' @examples #' set.seed(42) #' mat = matrix(runif(20), nrow = 5) -#' swm1 = matrix(runif(25), 5, 5) -#' swm2 = matrix(runif(25), 5, 5) +#' swm1 = apply(matrix(runif(25), 5, 5), 1, \(.x) .x / sum(.x)) +#' swm2 = apply(matrix(runif(25), 5, 5), 1, \(.x) .x / sum(.x)) #' coupling::metacoupling(mat, swm1, swm2) #' metacoupling = \(data, swm_peri = NULL, swm_tele = NULL, weight = NULL, From 1383ef9a9193c620bf3b0c7e593c8ca3eabd114e Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Thu, 18 Jun 2026 00:00:01 +0800 Subject: [PATCH 25/26] implement metacoupling wrappers --- man/metacoupling.Rd | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/man/metacoupling.Rd b/man/metacoupling.Rd index 1330fb5..f318dee 100644 --- a/man/metacoupling.Rd +++ b/man/metacoupling.Rd @@ -55,8 +55,8 @@ typically symmetric. \examples{ set.seed(42) mat = matrix(runif(20), nrow = 5) -swm1 = matrix(runif(25), 5, 5) -swm2 = matrix(runif(25), 5, 5) +swm1 = apply(matrix(runif(25), 5, 5), 1, \(.x) .x / sum(.x)) +swm2 = apply(matrix(runif(25), 5, 5), 1, \(.x) .x / sum(.x)) coupling::metacoupling(mat, swm1, swm2) } From 1553956cad8152add2654e52d1e1da7a2f83ea2b Mon Sep 17 00:00:00 2001 From: SpatLyu Date: Thu, 18 Jun 2026 00:00:01 +0800 Subject: [PATCH 26/26] implement metacoupling wrappers --- _pkgdown.yml | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/_pkgdown.yml b/_pkgdown.yml index b3d989a..36d1273 100644 --- a/_pkgdown.yml +++ b/_pkgdown.yml @@ -30,3 +30,7 @@ reference: - subtitle: coupling coordination degree contents: - ccd + +- subtitle: meta-coupling analysis + contents: + - metacoupling