scran_markers
Marker detection for single-cell data
Loading...
Searching...
No Matches
score_markers_pairwise.hpp
Go to the documentation of this file.
1#ifndef SCRAN_MARKERS_SCORE_MARKERS_PAIRWISE_HPP
2#define SCRAN_MARKERS_SCORE_MARKERS_PAIRWISE_HPP
3
4#include "cohens_d.hpp"
5#include "simple_diff.hpp"
6
7#include "scan_matrix.hpp"
8#include "average_group_stats.hpp"
9#include "PrecomputedPairwiseWeights.hpp"
10#include "create_combinations.hpp"
11
12#include <vector>
13
15#include "tatami/tatami.hpp"
16#include "tatami_stats/tatami_stats.hpp"
17
23namespace scran_markers {
24
76
81template<typename Stat_>
88 std::vector<Stat_*> mean;
89
95 std::vector<Stat_*> detected;
96
108 Stat_* cohens_d = NULL;
109
116 Stat_* auc = NULL;
117
124 Stat_* delta_mean = NULL;
125
132 Stat_* delta_detected = NULL;
133};
134
139template<typename Stat_>
145 std::vector<std::vector<Stat_> > mean;
146
151 std::vector<std::vector<Stat_> > detected;
152
159 std::vector<Stat_> cohens_d;
160
167 std::vector<Stat_> auc;
168
175 std::vector<Stat_> delta_mean;
176
183 std::vector<Stat_> delta_detected;
184};
185
189namespace internal {
190
191template<typename Index_, typename Stat_>
192void process_simple_pairwise_effects(
193 Index_ ngenes,
194 size_t ngroups,
195 size_t nblocks,
196 size_t ncombos,
197 std::vector<Stat_>& combo_means,
198 std::vector<Stat_>& combo_vars,
199 std::vector<Stat_>& combo_detected,
201 const std::vector<Stat_>& combo_weights,
202 double threshold,
203 int num_threads)
204{
205 std::vector<Stat_> total_weights_per_group;
206 const Stat_* total_weights_ptr = combo_weights.data();
207 if (nblocks > 1) {
208 total_weights_per_group = compute_total_weight_per_group(ngroups, nblocks, combo_weights.data());
209 total_weights_ptr = total_weights_per_group.data();
210 }
211 PrecomputedPairwiseWeights<Stat_> preweights(ngroups, nblocks, combo_weights.data());
212
213 tatami::parallelize([&](size_t, Index_ start, Index_ length) -> void {
214 size_t in_offset = ncombos * static_cast<size_t>(start);
215 const auto* tmp_means = combo_means.data() + in_offset;
216 const auto* tmp_variances = combo_vars.data() + in_offset;
217 const auto* tmp_detected = combo_detected.data() + in_offset;
218
219 size_t squared = ngroups * ngroups;
220 size_t out_offset = start * squared;
221 for (size_t gene = start, end = start + length; gene < end; ++gene, out_offset += squared) {
222 average_group_stats(gene, ngroups, nblocks, tmp_means, tmp_detected, combo_weights.data(), total_weights_ptr, output.mean, output.detected);
223
224 // Computing the effect sizes.
225 if (output.cohens_d != NULL) {
226 internal::compute_pairwise_cohens_d(tmp_means, tmp_variances, ngroups, nblocks, preweights, threshold, output.cohens_d + out_offset);
227 }
228
229 if (output.delta_detected != NULL) {
230 internal::compute_pairwise_simple_diff(tmp_detected, ngroups, nblocks, preweights, output.delta_detected + out_offset);
231 }
232
233 if (output.delta_mean != NULL) {
234 internal::compute_pairwise_simple_diff(tmp_means, ngroups, nblocks, preweights, output.delta_mean + out_offset);
235 }
236
237 tmp_means += ncombos;
238 tmp_variances += ncombos;
239 tmp_detected += ncombos;
240 }
241 }, ngenes, num_threads);
242}
243
244template<typename Stat_>
245ScoreMarkersPairwiseBuffers<Stat_> fill_pairwise_results(size_t ngenes, size_t ngroups, ScoreMarkersPairwiseResults<Stat_>& store, const ScoreMarkersPairwiseOptions& opt) {
246 ScoreMarkersPairwiseBuffers<Stat_> output;
247
248 internal::fill_average_results(ngenes, ngroups, store.mean, store.detected, output.mean, output.detected);
249
250 size_t num_effect_sizes = ngenes * ngroups * ngroups; // already size_t's, no need to cast.
251
252 if (opt.compute_cohens_d) {
253 store.cohens_d.resize(num_effect_sizes);
254 output.cohens_d = store.cohens_d.data();
255 }
256 if (opt.compute_auc) {
257 store.auc.resize(num_effect_sizes);
258 output.auc = store.auc.data();
259 }
260 if (opt.compute_delta_mean) {
261 store.delta_mean.resize(num_effect_sizes);
262 output.delta_mean = store.delta_mean.data();
263 }
264 if (opt.compute_delta_detected) {
265 store.delta_detected.resize(num_effect_sizes);
266 output.delta_detected = store.delta_detected.data();
267 }
268
269 return output;
270}
271
272}
341template<typename Value_, typename Index_, typename Group_, typename Stat_>
344 const Group_* group,
347{
348 Index_ NC = matrix.ncol();
349 auto group_sizes = tatami_stats::tabulate_groups(group, NC);
350
351 // In most cases this doesn't really matter, but we do it for consistency with the 1-block case,
352 // and to account for variable weighting where non-zero block sizes get zero weight.
353 auto group_weights = scran_blocks::compute_weights<Stat_>(group_sizes, options.block_weight_policy, options.variable_block_weight_parameters);
354
355 size_t ngroups = group_sizes.size();
356 size_t payload_size = static_cast<size_t>(matrix.nrow()) * ngroups; // cast to size_t to avoid overflow.
358
359 if (output.auc != NULL || matrix.prefer_rows()) {
360 internal::scan_matrix_by_row<true>(
361 matrix,
362 ngroups,
363 group,
364 1,
365 static_cast<int*>(NULL),
366 ngroups,
367 NULL,
371 output.auc,
374 options.threshold,
375 options.num_threads
376 );
377
378 } else {
379 internal::scan_matrix_by_column(
380 matrix,
381 ngroups,
382 group,
387 options.num_threads
388 );
389 }
390
391 internal::process_simple_pairwise_effects(
392 matrix.nrow(),
393 ngroups,
394 1,
395 ngroups,
399 output,
401 options.threshold,
402 options.num_threads);
403}
404
443template<typename Value_, typename Index_, typename Group_, typename Block_, typename Stat_>
446 const Group_* group,
447 const Block_* block,
450{
451 Index_ NC = matrix.ncol();
452 size_t ngroups = output.mean.size();
453 size_t nblocks = tatami_stats::total_groups(block, NC);
454
455 auto combinations = internal::create_combinations(ngroups, group, block, NC);
456 auto combo_sizes = internal::tabulate_combinations<Index_>(ngroups, nblocks, combinations);
457 size_t ncombos = combo_sizes.size();
458 auto combo_weights = scran_blocks::compute_weights<Stat_>(combo_sizes, options.block_weight_policy, options.variable_block_weight_parameters);
459
460 size_t payload_size = static_cast<size_t>(matrix.nrow()) * ncombos; // cast to size_t to avoid overflow.
462
463 if (output.auc != NULL || matrix.prefer_rows()) {
464 internal::scan_matrix_by_row<false>(
465 matrix,
466 ngroups,
467 group,
468 nblocks,
469 block,
470 ncombos,
471 combinations.data(),
475 output.auc,
478 options.threshold,
479 options.num_threads
480 );
481
482 } else {
483 internal::scan_matrix_by_column(
484 matrix,
485 ncombos,
486 combinations.data(),
491 options.num_threads
492 );
493 }
494
495 internal::process_simple_pairwise_effects(
496 matrix.nrow(),
497 ngroups,
498 nblocks,
499 ncombos,
503 output,
505 options.threshold,
506 options.num_threads);
507}
508
524template<typename Stat_ = double, typename Value_, typename Index_, typename Group_>
526 size_t ngroups = tatami_stats::total_groups(group, matrix.ncol());
528 auto buffers = internal::fill_pairwise_results(matrix.nrow(), ngroups, res, options);
530 return res;
531}
532
551template<typename Stat_ = double, typename Value_, typename Index_, typename Group_, typename Block_>
559
560}
561
562#endif
Marker detection for single-cell data.
Definition score_markers_pairwise.hpp:23
void score_markers_pairwise(const tatami::Matrix< Value_, Index_ > &matrix, const Group_ *group, const ScoreMarkersPairwiseOptions &options, const ScoreMarkersPairwiseBuffers< Stat_ > &output)
Definition score_markers_pairwise.hpp:342
void score_markers_pairwise_blocked(const tatami::Matrix< Value_, Index_ > &matrix, const Group_ *group, const Block_ *block, const ScoreMarkersPairwiseOptions &options, const ScoreMarkersPairwiseBuffers< Stat_ > &output)
Definition score_markers_pairwise.hpp:444
void parallelize(Function_ fun, Index_ tasks, int threads)
Buffers for score_markers_pairwise() and friends.
Definition score_markers_pairwise.hpp:82
Stat_ * delta_mean
Definition score_markers_pairwise.hpp:124
Stat_ * cohens_d
Definition score_markers_pairwise.hpp:108
std::vector< Stat_ * > detected
Definition score_markers_pairwise.hpp:95
Stat_ * delta_detected
Definition score_markers_pairwise.hpp:132
Stat_ * auc
Definition score_markers_pairwise.hpp:116
std::vector< Stat_ * > mean
Definition score_markers_pairwise.hpp:88
Options for score_markers_pairwise() and friends.
Definition score_markers_pairwise.hpp:28
scran_blocks::WeightPolicy block_weight_policy
Definition score_markers_pairwise.hpp:68
int num_threads
Definition score_markers_pairwise.hpp:39
bool compute_delta_mean
Definition score_markers_pairwise.hpp:57
double threshold
Definition score_markers_pairwise.hpp:33
bool compute_cohens_d
Definition score_markers_pairwise.hpp:45
bool compute_auc
Definition score_markers_pairwise.hpp:51
bool compute_delta_detected
Definition score_markers_pairwise.hpp:63
scran_blocks::VariableWeightParameters variable_block_weight_parameters
Definition score_markers_pairwise.hpp:74
Results for score_markers_pairwise() and friends.
Definition score_markers_pairwise.hpp:140
std::vector< Stat_ > delta_mean
Definition score_markers_pairwise.hpp:175
std::vector< std::vector< Stat_ > > detected
Definition score_markers_pairwise.hpp:151
std::vector< Stat_ > cohens_d
Definition score_markers_pairwise.hpp:159
std::vector< Stat_ > auc
Definition score_markers_pairwise.hpp:167
std::vector< Stat_ > delta_detected
Definition score_markers_pairwise.hpp:183
std::vector< std::vector< Stat_ > > mean
Definition score_markers_pairwise.hpp:145