LoginSignup
3

More than 5 years have passed since last update.

RedSVD を用いて疎行列の特異値分解を求める

Posted at

RedSVD

RedSVD は修正BSDライセンスで公開されている行列分解ライブラリです.

C++の線形代数ライブラリである Eigen を利用しており,簡単に利用することができます.
Eigenは疎行列を扱うためのクラスを提供していますが(SparseMatrix),特異値分解は密行列のみにしか対応しておりません.
RedSVDは SparseMatrix の特異値分解を計算することができるため,非常に便利です.

Wikiによるとheaderだけの RedSVDも公開されているため,今回はこちらのheaderを用いて実装してみたいと思います.

headerは /usr/local/include などに配置すれば利用することができます.

また,RedSVDや特異値分解の詳細は以下の記事を御覧ください.
DO++ 「行列分解ライブラリredsvdを公開しました」
http://hillbig.cocolog-nifty.com/do/2010/08/redsvd-aa59.html

実装例

実装例として,6×5の疎行列 $M$ の特異値分解を求めるプログラムを実装しました.

具体的には,$M = USV^T$ を満たす行列 $U$,$S$,$V$を求めます.

行列 $M$ は以下の通りです.

1 0 1 0 0
0 1 0 0 0
1 1 0 0 0
1 0 0 0 1
1 0 0 0 1
0 0 0 1 0

Eigenでは,Eigen::Triplet を用いて疎行列 SparseMatrix を作成することが一般的です.
vector<Eigen::Triplet<value_type>> に非零要素の位置と値を挿入し,SparseMatrix に対して setFromTriplets() を実行します.

このようにして疎行列を作成すれば,あとは RedSVD のコンストラクタを用いるだけで行列 $U$,$S$,$V$ を求めることができます.

RedSVD::RedSVD<Eigen::SparseMatrix<double>> svd(M);

コンストラクタ内の compute() によって既に特異値分解が計算されているので,あとは matrixU() といった関数で行列を取り出すだけです.

行列 $U$ は matrixU(),行列 $S$ は singularValues().asDiagonal(),行列 $V$ は matrixV() を実行します.

ここで注意なのですが,singularValues() は特異値が格納されたベクトルが得られるため,asDiagonal() 関数を用いて対角行列に変換する必要があります.

以上を踏まえ,疎行列 $M$ の特異値分解を求めるプログラムを以下に示します.

#include <iostream>
#include <vector>
#include <Eigen/Sparse>
#include <Eigen/SVD>
#include <RedSVD/RedSVD-h>

using namespace std;
using T = Eigen::Triplet<double>;

int main()
{
    Eigen::SparseMatrix<double> M(6, 5);
    vector<T> t;
    t.push_back(T(0, 0, 1));
    t.push_back(T(0, 2, 1));
    t.push_back(T(1, 1, 1));
    t.push_back(T(2, 0, 1));
    t.push_back(T(2, 1, 1));
    t.push_back(T(3, 0, 1));
    t.push_back(T(3, 4, 1));
    t.push_back(T(4, 0, 1));
    t.push_back(T(4, 4, 1));
    t.push_back(T(5, 3, 1));

    M.setFromTriplets(begin(t), end(t));

    cout << M << endl;

    RedSVD::RedSVD<Eigen::SparseMatrix<double>> svd(M);
    Eigen::MatrixXd U = svd.matrixU();
    Eigen::MatrixXd S = svd.singularValues().asDiagonal();
    Eigen::MatrixXd V = svd.matrixV();
    cout << "U: " << endl << U << endl;
    cout << "S: " << endl << S << endl;
    cout << "V: " << endl << V << endl;
    cout << "M_approx: " << endl << U*S*V.transpose() << endl;
}

実行結果

redsvd.png

まとめ

RedSVD を用いることで,簡単に疎行列の特異値分解を計算することができました.

今回の例では小さな行列に対して処理を行いましたが,より大きな行列に対しても実行可能です.
少しでも実験や開発の参考になれば幸いです.

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
3