2011年7月31日日曜日

正則な確率行列の逆行列

確率行列とは,行和が 1 な非負の正方行列のことです.確率行列 X に対して,Xi j を状態 i から状態 j への遷移確率と見なせば,X はマルコフ連鎖の遷移確率と見なすこともできます.

今正則な確率行列 X が与えられたとします.X の逆行列 X-1 は行和が 1 です.証明は簡単で,全成分が 1 である列ベクトルを 1 と書くことにすると,X の行和が 1 であることから 1 = X1 で,この式の両辺に左から X-1 をかければ X-11 = 1 となります.

これだけでは X-1 が確率行列とは言えません.一般に成分が負になる場合があるからです.逆行列が確率行列になるための条件は何でしょうか? とりあえず考えつく方法として 3 つの解法があるようです.

一番簡単なのは,X と X-1 を遷移確率行列だと思うことです.X が非決定的な遷移確率行列 (i に対して複数の j で正な値を取ることがある) の場合,X で遷移したあと,X-1 で遷移すると,確率 1 でもとの場所に戻ってくる,などという奇妙なことは絶対に起こりません.つまり,逆行列が確率行列になるのは,X が決定的な遷移確率行列の場合だけです.行列の言葉で言い換えると X が置換行列ということになります.

正面から証明する方法もあります.X-1 が非負行列で,X のある行 i に対して 2 つの列 j, j' で正の値をとっていると仮定します.XX-1 は単位行列なので,i 以外のすべてのインデックス i' に対して,X の i 行目と X-1 の i' 列目は内積が 0 です.これらは共に非負ベクトルなので,X-1 の (j, i'), (j', i') 成分はともに 0 となります.つまり X-1 の j, j' 行目は第 i 成分を除いて 0 となり,この 2 行は一次従属になります.これは X-1 が正則であることに反します.よって X はどの行も 1 つの成分を除いてゼロであり,置換行列になります.

最後に,固有値からも証明できるようです.確率行列はスペクトル半径が 1 です.X-1 は X の固有値の逆数を固有値に持つので,これらが共に確率行列となるためには,両方ともすべての固有値が長さ 1 である必要があります.つまり X は等長変換を表す行列 (直交行列) ということです.第 i 成分だけが 1 の行ベクトル v を考えると,vX のユークリッドノルムが 1 である必要があります.つまり X の第 i 行はユークリッドノルムが 1 となりますが,行和が 1 な非負ベクトルでユークリッドノルムも 1 となるのは,ある一つの成分だけが 1 というベクトル,しかありません.よって X は置換行列となります.


とまあ,きっかけがあってそんなことをちょっと考えていました.

2011年7月2日土曜日

eLog - C++ で簡易ログと簡単時間計測

C++ で使える簡単なロガーを書きました.ロガーと書いてしまいましたが,そんな大げさなものではなくて,ほとんど std::cerr の代わりみたいなものです.見た目はほとんど google-glog です.

beam2d/elog - GitHub

ヘッダーだけなので,すぐ使えます./usr/local/include へのインストールは ./waf install でできます.もちろん,自分のコードツリーにコピーしてきてもいいです.MIT License です.

使い方.

#include <elog/elog.h>

class ModuleA {};
class ModuleB {};

class MyLogger : public LOG::Logger { ... };

int main() {
  LOG(INFO) << "基本は google-glog と同じ";
  LOG(ERROR) << "LogLevel (severity) は INFO, WARN, ERROR, FATAL";

  try {
    LOG(FATAL) << "FATAL は";
  } catch (const LOG::FatalLogError&) {
    LOG(ERROR) << "例外を投げる";
  }

  LOG() << "INFO は省略できる";

  // ロガーのレベルを変更できる
  // LOG は括弧を付けなければ名前空間
  LOG::SetDefaultLoggerLevel(LOG::ERROR);

  LOG(INFO) << "出力されない";
  LOG(WARN) << "出力されない";

  LOG::SetDefaultLoggerLevel(LOG::INFO);

  LOG(ModuleA, 0) << "型ごとに verbosity をつけてログを吐ける";

  // 型ごとに異なる verbosity を設定できる (初期値は 0)
  LOG::SetDefaultLoggerVerbosity<ModuleA>(2);  // verbosity 2 以下を有効に
  LOG(ModuleA, 2) << "出力される";
  LOG(ModuleA, 3) << "出力されない";
  LOG(ModuleB, 1) << "出力されない";

  unsigned char uc = 65;
  LOG() << "(un)signed char は整数として出力される: " << uc;
  LOG() << "char は普通のストリームと同様: " << 'A';

  // アサートもある
  CHECK(true) << "出力されない";
  try {
    CHECK(false) << "出力され,そして";
  } catch (const LOG::CheckError&) {
    LOG(ERROR) << "例外を投げる";
  }

  // デフォルトのロガーは std::clog に吐くが,他のストリームに変えるには:
  std::ostream& ostream = std::cout;
  LOG::StreamLogger cout_logger(ostream);
  LOG::SetLogger(cout_logger);

  // この場合 SetDefaultLoggerLevel 等は使わない
  cout_logger.set_level(LOG::WARN);
  cout_logger.SetTypeVerbosity<ModuleA>(1);

  // ロガーの実装は自分で作って差し替えられる (Logger の定義は logger.h)
  MyLogger my_logger;
  LOG::SetLogger(my_logger);
}

ファイル名と行数も出力します.そのうち日付時刻とかも入れたいです.

あと簡単なベンチマーク機能があります.アルゴリズムをパパっと実装して速度を見るのに便利です.こっちは iwiwi くんの benchmark を参考にしています.

#include <elog/benchmark.h>

// BENCHMARK のログ出力はすべて INFO レベルです
int main() {
  // タイトルを出してから経過時間を計測して,最後に出力
  BENCHMARK(my_bench) {
    for (int i = 0; i < 100000; ++i);
    LOG() << "途中でも出力できる: " << my_bench.GetTime() << " sec";
    for (int i = 0; i < 100000; ++i);
  }

  // 複数の計測結果をまとめる
  LOG::BenchmarkSuite suite("title");

  BENCHMARK(suite, case_1) {
    for (int i = 0; i < 100000; ++i);
    LOG() << "こっちも同じく: " << case_1.GetTime() << " sec";
    for (int i = 0; i < 100000; ++i);
  }
  BENCHMARK(suite, case_2, "開始時にメッセージを出力できる") {
    for (int i = 0; i < 200000; ++i);
  }

  // 計測結果を表にして出力
  suite.LogChart();
}
とすると出力は
my_bench: ...
[INFO] bench.cc(8): 途中でも出力できる: 0.000623941 sec
my_bench: 0.00120592 sec
case_1: ...
[INFO] bench.cc(17): こっちも同じく: 0.000519037 sec
case_1: 0.00106406 sec
case_2: 開始時にメッセージを出力できる...
case_2: 0.00103211 sec
title   | time (sec)
--------+-----------
case_1  |    0.00106
case_2  |    0.00103
--------+-----------
(other) |    0.00006
total   |    0.00216

ログの構文は,周りの邪魔にならない簡素なものがいいと思っていて,google-glog はその意味で好きですが,マクロです.マクロは名前空間代わりのプレフィックスを付けた方が親切ですが,ログの構文は簡素な方がいいので,上のようになっています.その代わり,そういうグローバルスコープの単純な名前は極力減らす方針で書いています.

あと google-glog はリンクがいるし,最初に Init なんたらを呼ばなきゃいけないし,gflags との関係がややこしいので,何も考えずに使えるものが欲しかったということもあります.

型ごとに verbosity を設定するのはほぼ思いつきなので,便利なのかどうかわからないです.

というわけで,大したものではないですが,ぜひ.