computing column sums of matrix vector > w

Question

In a previous post column vector with row means -- with std::accumulate? I asked if it was possible, using STL functionality, to compute row means of a matrix

vector< vector<double> > data ( rows, vector<double> ( columns ) );

The top answer by @benjaminlindley is not only just what I was looking for, it is a thing of beauty. Forever hopeful I thought it would be as easy to compute column means, so an STL equivalent of

vector<double> colmeans( data[0].size() );
    for ( int i=0; i<data.size(); i++ )
        for ( int j=0; j<data[i].size(); j++ )            
            colmeans[j] += data[i][j]/data.size();

where the mean is not computed inside each vector<double>, but across the same index in all the vectors:

colmeans[0]       == ( data[0][0] + data[1][0] + ... data[rows][0] ) / rows
colmeans[1]       == ( data[0][1] + data[1][1] + ... data[rows][1] ) / rows
colmeans[2]       == ( data[0][2] + data[1][2] + ... data[rows][2] ) / rows
...
colmeans[columns] == ( data[0]   [columns] + 
                       data[1]   [columns] + 
                       ... 
                       data[rows][columns] ) / rows

It turns out to be quite different -- accumulate does not want to work on vectors of vectors. Is it somehow possible using accumulate with the [] operator? I cannot even come up with an intermediate form (to get rid of eather the for i or for j loop) which does not seem right.

Something with accumulate and the [] operator? Or bind?

Answer 1

Here's something I came up with, using for_each and transform:

std::vector<std::vector<double>> data { {1,2,3}, {1,2,3}, {1,2,3} };

std::vector<double> colsums( data[0].size() ); // initialize the size
                                                // to number of columns

std::for_each(data.begin(), data.end(),

    [&](const std::vector<double>& row)
    {
        // Use transform overload that takes two input ranges.
        // Note that colsums is the second input range as well as the output range.
        // We take each element of the row and add it to the corresponding
        // element of colsums vector:
        std::transform(row.begin(), row.end(), colsums.begin(), colsums.begin(),
                       [](double d1, double d2) { return d1 + d2; });
    });

std::cout << "Column means: ";
std::transform(
    colsums.begin(), colsums.end(),
    std::ostream_iterator<double>(std::cout, " "),
    [&data](double d) { return d / data.size(); });

LWS Demo

Answer 2

First let me state that you really should not nest std::vectors. Apart from that I got some solution which is of course longer than your initial code but that could save in the long run:

#include <vector>
#include <boost/iterator/iterator_adaptor.hpp>
#include <boost/iterator/counting_iterator.hpp>

typedef std::vector<std::vector<double> > Data;

struct ColumnElement : boost::iterator_adaptor<ColumnElement,
                                                Data::const_iterator,
                                                const double> {
        int col;

        ColumnElement(int col, const Data::const_iterator &iter)
        : iterator_adaptor(iter), col(col)
        {}
        const double& dereference()const { return (*base())[col]; }
};

struct Column {
        int col;
        const Data *data;

        Column(int col, const Data *data) : col(col), data(data) {}
        ColumnElement begin()const { return ColumnElement(col, data->begin()); }
        ColumnElement end()const { return ColumnElement(col, data->end()); }
        int size()const { return std::distance(begin(), end()); }
};

struct Columns : boost::iterator_adaptor<Columns, boost::counting_iterator<int>,
                                        Column, boost::use_default, Column> {
        const Data *data;

        Columns(int col, const Data *data): iterator_adaptor(col), data(data) {}

        Column dereference()const { return Column(*base(), data); }
};

Columns columnsBegin(const Data &data) { return Columns(0, &data); }
Columns columnsEnd(const Data &data) {
        return Columns(data.empty() ? 0 : data.front().size(), &data);
}

This could be used in short:

double Mean(const Column &d) {
        return std::accumulate(d.begin(), d.end(), 0.0) / d.size();
}

int main() {
        Data data = {   {1, 2, 3},
                        {2, 2, 2},
                        {9, 8, 7}};
        std::vector<double> colMeans(data[0].size());
        std::transform(columnsBegin(data), columnsEnd(data), 
                       colMeans.begin(), Mean);
        std::copy(colMeans.begin(), colMeans.end(),
                  std::ostream_iterator<double>(std::cout, ","));
        std::cout << "\n";
}

I employed some boost functionality to shorten it, but it could be done without boost (however much longer).

The idea was to create an iterator over all columns (called Columns just for shortness) and an iterator that iterates over all elements of one column (ColumnElement, also shortened, should be better named ColumnElementIterator) and Column that represents the range of all Elements of one column.