Added variant of combine_factors that reports unused combinations.

LTLA · LTLA · commit 9057609d7e05 · 2024-07-19T00:38:07.000-07:00
diff --git a/include/scran_aggregate/combine_factors.hpp b/include/scran_aggregate/combine_factors.hpp
@@ -93,6 +93,86 @@ std::vector<std::vector<Factor_> > combine_factors(size_t n, const std::vector<c
     return output;
 }
 
+/**
+ * This function is a variation of `combine_factors()` that considers unobserved combinations of factor levels.
+ *
+ * @tparam Factor_ Factor type.
+ * Any type may be used here as long as it is comparable.
+ * @tparam Number_ Integer type for the number of levels in each factor.
+ * @tparam Combined_ Integer type for the combined factor.
+ * This should be large enough to hold the number of unique (possibly unused) combinations.
+ *
+ * @param n Number of observations (i.e., cells).
+ * @param[in] factors Vector of pairs, each of which corresponds to a factor.
+ * The first element of the pair is a pointer to an array of length `n`, containing the factor level for each observation.
+ * The second element is the total number of levels for this factor, which may be greater than the largeset observed level.
+ * @param[out] combined Pointer to an array of length `n` in which the combined factor is to be stored.
+ * On output, each entry determines the corresponding observation's combination of levels by indexing into the inner vectors of the returned object;
+ * see the argument of the same name in `combine_factors()` for more details.
+ *
+ * @return 
+ * Vector of vectors containing each unique combinations of factor levels.
+ * This has the same structure as the output of `combine_factors()`,
+ * with the only difference being that unobserved combinations are also reported.
+ */
+template<typename Factor_, typename Number_, typename Combined_>
+std::vector<std::vector<Factor_> > combine_factors_unused(size_t n, const std::vector<std::pair<const Factor_*, Number_> >& factors, Combined_* combined) {
+    size_t nfac = factors.size();
+    std::vector<std::vector<Factor_> > output(nfac);
+
+    if (nfac > 1) {
+        // We iterate from back to front, where the first factor is the slowest changing.
+        std::copy_n(factors[nfac - 1].first, n, combined); 
+        Combined_ mult = factors[nfac - 1].second;
+        for (size_t f = nfac - 1; f > 0; --f) {
+            const auto& finfo = factors[f - 1];
+            auto ff = finfo.first;
+            for (size_t i = 0; i < n; ++i) {
+                combined[i] += mult * ff[i];
+            }
+            mult *= finfo.second;
+        }
+
+        auto ncombos = mult;
+        Combined_ outer_repeats = mult;
+        Combined_ inner_repeats = 1;
+        for (size_t f = nfac; f > 0; --f) {
+            auto& out = output[f - 1];
+            out.reserve(ncombos);
+
+            const auto& finfo = factors[f - 1];
+            size_t initial_size = inner_repeats * finfo.second;
+            out.resize(initial_size);
+
+            if (inner_repeats == 1) {
+                std::iota(out.begin(), out.end(), static_cast<Combined_>(0));
+            } else {
+                auto oIt = out.begin();
+                for (Number_ l = 0; l < finfo.second; ++l) {
+                    std::fill_n(oIt, inner_repeats, l);
+                    oIt += inner_repeats;
+                }
+            }
+            inner_repeats = initial_size;
+
+            outer_repeats /= finfo.second;
+            for (Combined_ r = 1; r < outer_repeats; ++r) {
+                out.insert(out.end(), out.begin(), out.begin() + initial_size);
+            }
+        }
+
+    } else if (nfac == 1) {
+        output[0].resize(factors[0].second);
+        std::iota(output[0].begin(), output[0].end(), static_cast<Combined_>(0));
+        std::copy_n(factors[0].first, n, combined);
+
+    } else {
+        std::fill_n(combined, n, 0);
+    }
+
+    return output;
+}
+
 }
 
 #endif
diff --git a/tests/src/combine_factors.cpp b/tests/src/combine_factors.cpp
@@ -45,6 +45,13 @@ TEST(CombineFactors, Simple) {
         std::vector<int> levels { 1, 3, 5, 7, 9 };
         EXPECT_EQ(combined.first[0], levels);
     }
+
+    // Nothing at all.
+    {
+        auto combined = test_combine_factors(10, std::vector<const int*>{});
+        EXPECT_EQ(combined.second, std::vector<int>(10));
+        EXPECT_TRUE(combined.first.empty());
+    }
 }
 
 TEST(CombineFactors, Multiple) {
@@ -101,3 +108,93 @@ TEST(CombineFactors, Multiple) {
         EXPECT_EQ(factor2, combined.first[1]);
     }
 }
+
+template<typename Factor_, typename Number_>
+std::pair<std::vector<std::vector<Factor_> >, std::vector<int> > test_combine_factors_unused(size_t n, const std::vector<std::pair<const Factor_*, Number_> >& factors) {
+    std::vector<int> combined(n);
+    auto levels = scran_aggregate::combine_factors_unused(n, factors, combined.data());
+    return std::make_pair(std::move(levels), std::move(combined));
+}
+
+TEST(CombineFactorsUnused, Basic) {
+    std::vector<int> stuff{ 1, 3, 5, 3, 1 };
+    auto combined = test_combine_factors_unused(stuff.size(), std::vector<std::pair<const int*, int> >{ { stuff.data(), 7 } });
+    EXPECT_EQ(combined.second, stuff);
+    std::vector<int> levels{ 0, 1, 2, 3, 4, 5, 6 };
+    EXPECT_EQ(combined.first[0], levels);
+
+    auto combined2 = test_combine_factors_unused(10, std::vector<std::pair<const int*, int> >{});
+    EXPECT_EQ(combined2.second, std::vector<int>(10));
+    EXPECT_TRUE(combined2.first.empty());
+}
+
+TEST(CombineFactorsUnused, Multiple) {
+    std::vector<int> stuff1{ 0, 0, 1, 1, 1, 2, 2, 2, 2 };
+    std::vector<int> stuff2{ 0, 1, 2, 0, 1, 2, 0, 1, 2 };
+
+    {
+        auto combined = test_combine_factors_unused(stuff1.size(), 
+            std::vector<std::pair<const int*, int> >{ 
+                { stuff1.data(), 3 },
+                { stuff2.data(), 3 }
+            }
+        );
+
+        std::vector<int> expected { 0, 1, 5, 3, 4, 8, 6, 7, 8 };
+        EXPECT_EQ(combined.second, expected);
+
+        EXPECT_EQ(combined.first.size(), 2);
+        std::vector<int> levels1 { 0, 0, 0, 1, 1, 1, 2, 2, 2 };
+        std::vector<int> levels2 { 0, 1, 2, 0, 1, 2, 0, 1, 2 }; // (0, 2) is included even though it is not observed.
+        EXPECT_EQ(combined.first[0], levels1);
+        EXPECT_EQ(combined.first[1], levels2);
+    }
+
+    {
+        auto combined = test_combine_factors_unused(stuff1.size(), 
+            std::vector<std::pair<const int*, int> >{ 
+                { stuff1.data(), 4 },
+                { stuff2.data(), 3 }
+            }
+        );
+
+        std::vector<int> expected { 0, 1, 5, 3, 4, 8, 6, 7, 8 };
+        EXPECT_EQ(combined.second, expected);
+
+        EXPECT_EQ(combined.first.size(), 2);
+        std::vector<int> levels1 { 0, 0, 0, 1, 1, 1, 2, 2, 2, 3, 3, 3 };
+        std::vector<int> levels2 { 0, 1, 2, 0, 1, 2, 0, 1, 2, 0, 1, 2 }; // (0, 2) and (3, *) are included even though they are not observed.
+        EXPECT_EQ(combined.first[0], levels1);
+        EXPECT_EQ(combined.first[1], levels2);
+    }
+
+    // Multiple things at play here.
+    {
+        std::vector<int> mock{ 1 };
+        auto combined = test_combine_factors_unused(1,
+            std::vector<std::pair<const int*, int> >{ 
+                { mock.data(), 2 },
+                { mock.data(), 3 },
+                { mock.data(), 4 }
+            }
+        );
+
+        EXPECT_EQ(combined.second[0], 17); // i.e., 1*3*4 + 1*4 + 1.
+
+        auto create_mock_sequence = [](size_t nlevels, size_t inner, size_t outer) -> auto { 
+            std::vector<int> output(nlevels * inner * outer);
+            auto oIt = output.begin();
+            for (size_t o = 0; o < outer; ++o) {
+                for (size_t l = 0; l < nlevels; ++l) {
+                    std::fill_n(oIt, inner, l);
+                    oIt += inner;
+                }
+            }
+            return output;
+        };
+
+        EXPECT_EQ(combined.first[0], create_mock_sequence(2, 12, 1));
+        EXPECT_EQ(combined.first[1], create_mock_sequence(3, 4, 2));
+        EXPECT_EQ(combined.first[2], create_mock_sequence(4, 1, 6));
+    }
+}
