numpy そっくりなヘッダオンリー C++ ライブラリ NumCpp

はじめに

先日 Twitter でこんなツイートを見ました。

Pythonのいいところって大体ほとんどのライブラリがnumpyのオブジェクト使ってるところな気がする。C++もそういうのあればいいのに、あるかもしらんけど

— 🔰😢神😢🔰 (@AllisTheWill) September 3, 2019

確かに numpy はよく出来てて、Python のシンタックスにマッチしたとても便利なライブラリですね。このツイートのレスに以前見つけて遊んだ事のある NumCpp を「numpy そっくりさん」という内容で言及しました。

一応 numpy そっくりさんとして NumCpp もあります。 https://t.co/vamY0Q9qiy

— mattn＠有益情報 (@mattn_jp) September 3, 2019

ただ NumCpp がどんな使い勝手なのか知らない人も多いと思いますので、本記事では NumCpp を使ったサンプルを紹介したいと思います。

基本

numpy の np.array 相当は nc::NdArray のコンストラクタに該当します。NumCpp ではほぼ全て NdArray で処理を行います。NdArray はN次元の配列になっており、1次元であっても 1xN な配列として扱われます。

nc::NdArray<int> a = { {1, 2}, {3, 4}, {5, 6} };
std::cout << a << std::endl;

標準出力すると numpy 同様に表示されます。

[[1.000000, 2.000000, ]
[3.000000, 4.000000, ]
[5.000000, 6.000000, ]]

reshape

numpy 同様、reshape 時に内部でメモリが再確保されている訳でなく、要素サイズと要素数が変更されているだけなので気軽に形状を変更できます。

a.reshape(2, 3);

[[1.000000, 2.000000, 3.000000, ]
[4.000000, 5.000000, 6.000000, ]]

演算子

また演算子もオーバロードされていて、簡単に行列演算も行えます。

nc::NdArray<float> a = {{1, 2}, {3, 4}};
nc::NdArray<float> b = {{2, 3}, {4, 5}};
std::cout << (a + b) * 2 << std::endl;

[[6.000000, 10.000000, ]
[14.000000, 18.000000, ]]

乱数

乱数も簡単に生成できます。

auto a = nc::random::randN<double>({3, 4});
std::cout << a << std::endl;

[[-0.150035, -0.444410, -1.678120, 2.476207, ]
[-0.528850, 0.659141, -0.538446, 1.673758, ]
[-0.111288, 2.069479, 1.395632, -1.650037, ]]

注意点

一点、気を付けるべき点としては for で回すと次元単位ではなく要素単位になってしまう点です。

nc::NdArray<float> a = {{1, 2}, {3, 4}};
for(auto v : a) {
  std::cout << v << std::endl;
}

1
2
3
4

これは添え字を指定した場合も同じです。

nc::NdArray<float> a = {{1, 2}, {3, 4}};
a[2] = 5;
std::cout << a << std::endl;

[[1.000000, 2.000000, ]
[5.000000, 4.000000, ]]

これを次元単位にするには以下の様に書きます。

nc::NdArray<float> a = {{1, 2}, {3, 4}};
for(int i = 0; i < a.shape().rows; i++) {
  std::cout << a.row(i) << std::endl;
}

[[1.000000, 2.000000, ]]
[[3.000000, 4.000000, ]]

その他、numpy にある norm や dot も見た目そのままに使えます。詳しくは NumCpp の README を参照下さい。

ちょっとしたサンプル

最後に NumCpp を使って、あやめのCSVデータをロジスティック回帰で推論するサンプルを載せておきます。

#include <NumCpp.hpp>
#include <vector>
#include <map>
#include <string>
#include <iostream>

static float
softmax(nc::NdArray<float>& w, nc::NdArray<float>& x) {
  auto v = nc::dot<float>(w, x).item();
  return 1.0 / (1.0 + std::exp(-v));
}

static float
predict(nc::NdArray<float>& w, nc::NdArray<float>& x) {
  return softmax(w, x);
}

static nc::NdArray<float>
logistic_regression(nc::NdArray<float>& X, nc::NdArray<float>& y, float rate, int ntrains) {
  auto w = nc::random::randN<float>({1, X.shape().cols});

  for (int n = 0; n < ntrains; n++) {
    for (int i = 0; i < X.shape().rows; i++) {
      auto x = X.row(i);
      auto t = nc::NdArray<float>(x.copy());
      auto pred = softmax(t, w);
      auto perr = y.at(0, i) - pred;
      auto scale = rate * perr * pred * (1 - pred);
      auto dx = x.copy();
      dx += x;
      dx *= scale;
      w += dx;
    }
  }

  return w;
}

static std::vector<std::string>
split(std::string& fname, char delimiter) {
  std::istringstream f(fname);
  std::string field;
  std::vector<std::string> result;
  while (getline(f, field, delimiter)) {
      result.push_back(field);
  }
  return result;
}

int main() {
  std::ifstream ifs("iris.csv");

  std::string line;

  // skip header
  std::getline(ifs, line);

  std::vector<float> rows;
  std::vector<std::string> names;
  while (std::getline(ifs, line)) {
    // sepal length, sepal width, petal length, petal width, name
    auto cells = split(line, ',');
    rows.push_back(std::stof(cells.at(0)));
    rows.push_back(std::stof(cells.at(1)));
    rows.push_back(std::stof(cells.at(2)));
    rows.push_back(std::stof(cells.at(3)));
    names.push_back(cells.at(4));
  }
  // make vector 4 dimentioned
  nc::NdArray<float> X(rows);
  X.reshape(rows.size()/4, 4);

  // make onehot values of names
  std::map<std::string, int> labels;
  std::for_each(names.begin(), names.end(), [&](decltype(names)::value_type x) {
    if (labels.count(x) == 0) labels[x] = labels.size();
  });
  std::vector<float> n;
  std::for_each(names.begin(), names.end(), [&](decltype(names)::value_type x) {
    if (labels.count(x) > 0) n.push_back((float)labels[x]);
  });
  nc::NdArray<float> y(n);
  y /= labels.size();

  names.clear();
  for(auto k : labels) {
    names.push_back(k.first);
  }

  // make factor from input values
  auto w = logistic_regression(X, y, 0.1, 300);

  // predict samples
  for (int i = 0; i < X.shape().rows; i++) {
    auto x = X.row(i);
    auto n = (int) (predict(w, x) * (float) labels.size() + 0.5);
    std::cout << names[n] << std::endl;
  }

  return 0;
}

結論

NumCpp 便利。ただ blas 等で高速化されている訳ではないので、カリッカリにチューンされた環境で役立つかどうかは試してみないと分かりません。ご注意下さい。