compute.hpp

Go to the documentation of this file.

#ifndef GSDECON_COMPUTE_HPP
#define GSDECON_COMPUTE_HPP
 
#include <algorithm>
#include <numeric>
 
#include "tatami/tatami.hpp"
#include "irlba/irlba.hpp"
#include "scran_pca/scran_pca.hpp"
 
#include "Options.hpp"
#include "Results.hpp"
#include "utils.hpp"
 
namespace gsdecon {
 
template<typename Value_, typename Index_, typename Float_>
void compute(const tatami::Matrix<Value_, Index_>& matrix, const Options& options, const Buffers<Float_>& output) {
    if (internal::check_edge_cases(matrix, options.rank, output)) {
        return;
    }
 
    scran_pca::SimplePcaOptions sopt;
    sopt.number = options.rank;
    sopt.scale = options.scale;
    sopt.realize_matrix = options.realize_matrix;
    sopt.num_threads = options.num_threads;
    sopt.irlba_options = options.irlba_options;
    auto res = scran_pca::simple_pca(matrix, sopt);
 
    double shift = std::accumulate(res.center.begin(), res.center.end(), 0.0) / matrix.nrow();
    std::fill_n(output.scores, matrix.ncol(), shift);
    internal::process_output(res.rotation, res.components, options.scale, res.scale, output);
}
 
template<typename Float_ = double, typename Value_, typename Index_>
Results<Float_> compute(const tatami::Matrix<Value_, Index_>& matrix, const Options& options) {
    Results<Float_> output;
    output.weights.resize(matrix.nrow());
    output.scores.resize(matrix.ncol());
 
    Buffers<Float_> buffers;
    buffers.weights = output.weights.data();
    buffers.scores = output.scores.data();
 
    compute(matrix, options, buffers);
    return output;
}
 
}
 
#endif