👍📕😓🎉😆⚽✌
昨今、世界の砂漠化と並び、世界の絵文字(Emoji)化が進んでおり、社会問題として取り上げられています。
- 英オクスフォード辞書の「今年の言葉」のEmoji化
- プログラミング言語のEmoji化
- 画像のEmoji化
- チャットのEmoji化
- hashtagのEmoji化
- ニューヨーカー紙の表紙のEmoji化
- など
自然言語もEmoji化が進んでおり、今後数年以内に「日本語はEmojiにすべて置き換わる」とも言われております。
そこで、本ブログもこの流れを汲み、自然言語処理・日本語処理からEmoji処理・Emojineeringを中心とした技術へとシフトしていきたいと思います。
用いる言語を日本語からEmojiへｔ👍📕😓🎉😆⚽✌😙
ラティスのNbestを求める
はじめに
形態素解析とかの解析時に使うラティス(形態素候補をグラフにしたもの)のうち、1番ベストな解だけが欲しい場合が多いが、2番目以降の解も欲しい場合がある。
N番目までの解を効率よく求める方法があり、使いたいケースが出てきたのに書いてみる。
Nbestとは
- 解析したときに、スコアが良い順に、第一候補(1best)の解だけでなく、第二候補、第三候補・・・の解が考えられるとき、第N候補までのことをNbestという
- Nbest解
前向きDP後ろ向きA*アルゴリズム
- http://ci.nii.ac.jp/naid/110002725063
- ラティスにおいて、効率よくNbestを求める方法が提案されている
- 最初に、1bestを求める方法と同じようにスタートノードからあるノードまでの最小コストhをすべてのノードについて求めておく
- ゴールノードからスタートノードに向かって、今考えているノードから(最適とは限らない、これまで取ったパスによって変わる)ゴールまでのコストgを求めていく
- 前側のコスト合計hと後ろ側のコスト合計gを使うと、スタートから今考えているノードを通って(今まで通ってきた最適とは限らない)ゴールまでのパスの合計値fが求められる
- fを使って小さい順に考えるノードを増やしながら、スタートノードにたどり着くパスを見つけると、Nbestの解を順次求められる
書いてみたけどわかりにくい。
徳永さんの「日本語入力を支える技術」の中で丁寧に紹介されている。
コード
スタートノード(BOS)とゴールノード(EOS)が1つずつあり、必ずスタートノードからゴールノードにパスを持つようなDAGであれば大丈夫だと思うので、それで書いてみる。
(雑なので、パスのノードを別に持ったり、トポロジカル順序を求めてる無駄がある・・・)
#include <iostream> #include <vector> #include <queue> #include <algorithm> #include <map> #include <set> #define INF 99999999 struct Node { int id; int w; //id,w std::map<int,int> forward; std::map<int,int> backward; int h; }; struct Path { int id; int g; int f; std::vector<int> rev; //現在idからゴールノードまでの逆順パス(注意:無駄にメモリを使っている) }; bool operator>(const Path& a, const Path& b){ return a.f > b.f; } //スタートとゴールとなるノードがそれぞれ1つずつあり、ゴールノード以外で終わるようなノードが無いDAG class stDAG { int n; //ノード数 int s, t; //スタートノード番号、ゴールノード番号 std::vector<Node> nodes; //ノード情報 //トポロジカル順序を求めるため bool topological_sort(std::vector<int>& order){ std::vector<int> color(n); for(int i=0; i<n; i++){ if(!color[i] && !topological_visit(i, order, color)) return false; } std::reverse(order.begin(), order.end()); return true; } bool topological_visit(int v, std::vector<int>& order, std::vector<int>& color){ color[v] = 1; for(std::map<int,int>::iterator itr = nodes[v].forward.begin(); itr != nodes[v].forward.end(); ++itr){ if(color[itr->first] == 2) continue; if(color[itr->first] == 1) return false; if(!topological_visit(itr->first, order, color)) return false; } order.push_back(v); color[v] = 2; return true; } public: //ノード数、スタートノードID、ゴールノードID stDAG(int n, int s, int t): n(n), nodes(n), s(s), t(t) { for(int i=0; i<n; i++){ nodes[i].id = i; nodes[i].h = INF; } } //ノード情報の追加 void add_node(int id, int w){ nodes[id].w = w; } //エッジ情報の追加 void add_edge(int id, int next_id, int w){ nodes[id].forward[next_id] = w; nodes[next_id].backward[id] = w; } //前向きDP void forward_dp(){ std::vector<int> order; if(!topological_sort(order)){ std::cerr << "topological sort error" << std::endl; return; } if(nodes[order[0]].id != s || nodes[order[n-1]].id != t){ std::cerr << "s,t id error" << std::endl; return; } nodes[order[0]].h = 0; for(int i=0; i<n; i++){ int id = order[i]; //std::cout << id << "\t" << nodes[id].h << std::endl; for(std::map<int,int>::iterator itr = nodes[id].forward.begin(); itr != nodes[id].forward.end(); ++itr){ int next_id = itr->first; int next_w = nodes[next_id].w; int edge_w = itr->second; nodes[next_id].h = std::min(nodes[next_id].h, nodes[id].h + edge_w + next_w); } } } //後向きA* void backward_astar(int N){ int no = 1; std::priority_queue<Path, std::vector<Path>, std::greater<Path> > que; Path path; path.id = t; path.g = 0; path.f = nodes[t].h; path.rev.push_back(t); que.push(path); while(!que.empty()){ path = que.top(); que.pop(); int id = path.id; int bestf = path.f; if(id == s){ std::cout << no << "best: "; for(int i=path.rev.size()-1; i>=0; i--){ std::cout << path.rev[i] << " "; } std::cout << " => " << path.f + nodes[s].w << std::endl; if(no == N) return; no++; } for(std::map<int,int>::iterator itr = nodes[id].backward.begin(); itr != nodes[id].backward.end(); ++itr){ int gg = nodes[id].w + itr->second + path.g; int ff = gg + nodes[itr->first].h; Path new_path; new_path.id = itr->first; new_path.g = gg; new_path.f = ff; new_path.rev = path.rev; new_path.rev.push_back(itr->first); que.push(new_path); } } } }; int main(){ //graph1 stDAG dag(6, 0, 5); dag.add_node(0, 0); dag.add_node(1, 1); dag.add_node(2, 2); dag.add_node(3, 3); dag.add_node(4, 4); dag.add_node(5, 5); dag.add_edge(0, 1, 1); dag.add_edge(0, 2, 2); dag.add_edge(1, 3, 2); dag.add_edge(1, 4, 1); dag.add_edge(2, 3, 1); dag.add_edge(2, 4, 2); dag.add_edge(3, 5, 2); dag.add_edge(4, 5, 1); dag.forward_dp(); dag.backward_astar(4); //graph2 /* stDAG dag(7, 0, 6); dag.add_node(0, 3); dag.add_node(1, 1); dag.add_node(2, 2); dag.add_node(3, 1); dag.add_node(4, 2); dag.add_node(5, 3); dag.add_node(6, 2); dag.add_edge(0, 1, 2); dag.add_edge(0, 2, 1); dag.add_edge(1, 3, 3); dag.add_edge(2, 3, 1); dag.add_edge(3, 4, 1); dag.add_edge(3, 5, 2); dag.add_edge(4, 6, 4); dag.add_edge(5, 6, 1); dag.forward_dp(); dag.backward_astar(4); */ return 0; }
結果
グラフ1
$ ./a.out 1best: 0 1 4 5 => 13 2best: 0 1 3 5 => 14 3best: 0 2 3 5 => 15 4best: 0 2 4 5 => 16
グラフ2
$ ./a.out 1best: 0 2 3 5 6 => 16 2best: 0 2 3 4 6 => 17 3best: 0 1 3 5 6 => 18 4best: 0 1 3 4 6 => 19
とりあえず合ってそうなので、大丈夫そう。
pandas使ってwikipediaの表データを取得する
はじめに
特定ジャンルの用語などをまとめて取得するのに、wikipediaの「〜の一覧」が有用だったりする。
wikipedia:一覧の一覧
多くは、リスト形式で書かれていたりするが、中には表(テーブル)形式でまとめられているものもある。
いろんな取得方法が考えられるが、pythonのデータ解析支援ライブラリの「pandas」を使うと簡単に取得できたのでメモしておく。
注意
記事を大量にダウンロードするためにクローラを使わないで下さい。 強引なクローリングは、ウィキペディアが劇的に遅くなる原因となります。 ウィキペディアのデータベースから自動的にデータの収集がなされた場合、 システム管理者によってあなたのサイトからウィキペディアへのアクセスを 禁止する措置が取られることもあります。 またウィキメディア財団が法的措置を検討することもあります。
自己責任でお願いします。
表形式でまとめられている一覧ページの取得
「猫の品種の一覧」「犬の品種一覧」「ウサギの品種一覧」はそれぞれ写真付きで表形式一覧になっている。
DBpediaで見てみても、表形式のデータをきれいに取り出せてはいない模様。
以下で取得、ちょっとした置換処理、csv出力ができる。
import pandas # 猫の品種の一覧 url = 'https://ja.wikipedia.org/wiki/%E7%8C%AB%E3%81%AE%E5%93%81%E7%A8%AE%E3%81%AE%E4%B8%80%E8%A6%A7' df = pandas.io.html.read_html(url, encoding='utf-8', header=0) df[0]['種類'] = df[0]['種類'].replace({'（英語版）':'', '\(.+\)':'', '\[\d+\]':'', '^en\:':''}, regex=True) df[0].to_csv('neko.csv') # 犬の品種一覧 url = 'https://ja.wikipedia.org/wiki/%E7%8A%AC%E3%81%AE%E5%93%81%E7%A8%AE%E4%B8%80%E8%A6%A7' df = pandas.io.html.read_html(url, encoding='utf-8', header=[0,1]) df[0]['犬種'] = df[0]['犬種'].replace({'（英語版）':'', '\(.+\)':'', '\[\d+\]':'', '^en\:':''}, regex=True) df[0].to_csv('inu.csv') # ウサギの品種一覧 url = 'https://ja.wikipedia.org/wiki/%E3%82%A6%E3%82%B5%E3%82%AE%E3%81%AE%E5%93%81%E7%A8%AE%E4%B8%80%E8%A6%A7' df = pandas.io.html.read_html(url, encoding='utf-8', header=0) df[0]['和名'] = df[0]['和名'].replace({'（英語版）':'', '\(.+\)':'', '\[\d+\]':'', '^en\:':''}, regex=True) df[0]['英名'] = df[0]['英名'].replace({'\[\d+\]':'', '^en\:':''}, regex=True) df[0].to_csv('usagi.csv')
結果
$ python get_wiki.py
$ ls
get_wiki.py inu.csv neko.csv usagi.csv
$ wc -l *.csv
759 inu.csv
91 neko.csv
130 usagi.csv
980 total
$ head *.csv
==> inu.csv <==
,犬種,原産地,国際畜犬連盟[1],JKC登録,画像,Unnamed: 5_level_0
,分類,番号,Unnamed: 2_level_1,Unnamed: 3_level_1,Unnamed: 4_level_1,Unnamed: 5_level_1
,,,,,,
0,アーカンソー・ジャイアント・ブルドッグ,アメリカ合衆国,,,,
1,アーフェンピンシャー,ドイツ・フランス,2G-1,186,○,
2,アイスランド・シープドッグ,アイスランド,5G-3,289,,
3,アイディ,モロッコ,2G-2,247,,
4,アイリッシュ・ウォーター・スパニエル,アイルランド,8G-3,124,○,
5,アイリッシュ・ウルフハウンド,アイルランド,10G-2,160,○,
6,アイリッシュ・コリー,アイルランド,,,,
==> neko.csv <==
,種類,原産国,発生,身体のタイプ,毛足の長さ,毛色および模様,画像
0,アビシニアン,エチオピア,自然発生種,フォーリン,短毛,複色,
1,アメリカンカール,アメリカ合衆国,突然変異,セミフォーリン,短毛もしくは長毛,単色,
2,アメリカンショートヘア,アメリカ合衆国,自然発生種,セミコビー,短毛,All,
3,アメリカンボブテイル,アメリカ合衆国,自然発生種,ロング＆サブスタンシャル,短毛もしくは長毛,All,
4,アメリカンワイヤーヘア,アメリカ合衆国,突然変異,セミコビー,Rex,All but colorpoint,
5,エキゾチックショートヘア,アメリカ合衆国,交配種,コビー,Short,All,
6,エジプシャンマウ,エジプト,自然発生種,セミフォーリン,短毛,ぶち,
7,オシキャット,アメリカ合衆国,交配種,セミフォーリン,短毛,Spotted,
8,オホースアズーレス,アメリカ合衆国,突然変異,フォーリン,Short,All,
==> usagi.csv <==
,和名,英名,原産,体重(kg),毛,耳,色,ARBA,BRC,画像
0,,Alaska,ドイツ,3.2-4.1,短,立,黒,No,Yes,
1,,Altex,アメリカ,5.9,短,立,Pointed White,No,No,---
2,アメリカン,American Blue,アメリカ,4.1-5.4,短,立,青,Yes,No,
3,アメリカン,American Albino,アメリカ,4.1-5.4,短,立,白,Yes,No,
4,アメリカン・ファジー・ロップ,American Fuzzy Lop,アメリカ,1.4-1.8,長,垂,さまざま,Yes,No,
5,アメリカン・セーブル,American Sable,アメリカ,4.1-4.5,短,立,Sable,Yes,No,
6,,Argente Bleu,フランス,2.7,短,立,青,No,Yes,---
7,,Argente Brun,フランス,2.7,短,立,茶,No,Yes,---
8,,Argente Clair,フランス,2.7,短,立,Blue Silver,No,No,---合わせて1000種類程度の品種名が取得できた。
wikipediaに限った話ではないので、お手軽に表形式の情報抽出するならかなりよさそう。
参考
- http://qiita.com/kitsuyui/items/4906bb457af4d0e2d0a5
- Wes McKinney, Pythonによるデータ分析入門 NumPy、pandasを使ったデータ処理(日本語訳), O'Reillyジャパン
Minimal Acyclic Subsequential Transducerで遊ぶ
はじめに
https://pycon.jp/2015/ja/proposals/vote/11/
Pycon2015で発表された「Pythonで作って学ぶ形態素解析」で紹介されていた辞書データ構造の「Minimal Acyclic Subsequential Transducer」について、勉強のために書いてみた。
Minimal Acyclic Subsequential Transducerとは
- Finite State Transducerの一種
- Transducerにおいて、initial stateが一つで、同じ入力ラベルを共有する同じ状態からのの遷移が2つ以上なく、各最終状態での最終出力文字列が高々p個のとき、p-subsequentialで、pが整数ならfinitely subsequentialというらしい
- minimal(状態数が最少)、Acyclic(サイクルが無い)
- Cyril Allauzen, Mehryar Mohri, Finitely Subsequential Transducersとp-Subsequentiable Transducersあたりを読む
- 「Lucene/Solrで使われているFST」はこれの事
わかりやすい他の方の参考資料
構築方法
- Stoyan Mihov, Denis Maurel, Direct Construction of Minimal Acyclic Subsequential Transducers, 2001
OpenFstとか使うイメージでいたけど、上記で紹介されている通り、一時的な状態を作らずに、ソート済みの入力から直接FSTを構築する方法が提案されている。
詳しい手順も上記の資料(qiita)で紹介されいているので、省略。
コード
論文にできるだけ従って、入力・出力は文字列、出力は一つのみ、で実装。
問題がある部分は、状態の探索(member関数)で線形に等価な状態を探索している点や、状態の等価性判定(equal関数)が下記の方法でもよいのか怪しい点(Def.2-3？)、など。
他の実装は全然見てないので、間違ってたら後で修正する。
#include <iostream> #include <vector> #include <algorithm> #include <list> #include <cstdio> #include <map> #include <queue> //入力文字列の最大長 #define MAX_WORD_SIZE 100 struct State { State* next[0x100]; std::string output[0x100]; std::string state_out; bool final; int _num; State(){ clear(); } State(State* s){ _num = -1; state_out = s->state_out; final = s->final; for(size_t i=0; i<0x100; i++){ next[i] = s->next[i]; output[i] = s->output[i]; } } ~State(){} void clear(){ _num = -1; state_out = ""; final = false; for(size_t i=0; i<0x100; i++){ next[i] = (State*)0; output[i] = ""; } } //状態の等価チェック // 与えられた状態sに対して、遷移状態と出力記号などが完全に一致するかを確認 bool equal(State* s){ if(final && s->final){ if(state_out == s->state_out) return true; return false; } std::queue< std::pair< std::pair<State*,std::string>, std::pair<State*,std::string> > > que; for(size_t i=0; i<0x100; i++){ que.push(std::make_pair(std::make_pair(next[i],""), std::make_pair(s->next[i],""))); } while(!que.empty()){ std::pair<std::pair<State*,std::string>, std::pair<State*,std::string> > p = que.front(); que.pop(); if(p.first.second != p.second.second) return false; if(p.first.first == NULL && p.second.first == NULL) continue; if(p.first.first == NULL || p.second.first == NULL) return false; if(p.first.first->state_out != p.second.first->state_out) return false; for(size_t i=0; i<0x100; i++){ if(p.first.first->next[i] == NULL && p.second.first->next[i] == NULL) continue; que.push(std::make_pair(std::make_pair(p.first.first->next[i],p.first.second + (char)i), std::make_pair(p.second.first->next[i],p.second.second + (char)i))); } } return true; } }; struct Dictionary { Dictionary(){} ~Dictionary(){ for(std::list<State*>::iterator itr = states.begin(); itr != states.end(); ++itr){ delete *itr; } states.clear(); for(size_t i=0; i<MAX_WORD_SIZE; i++){ if(TempState[i]){ delete TempState[i]; TempState[i] = NULL; } } } //辞書の状態に存在する状態かどうかチェック // [注意] 線形に全部等価チェックしているので非常に重い State* member(State* s){ for(std::list<State*>::iterator itr = states.begin(); itr != states.end(); ++itr){ State* p = *itr; if(s->equal(p)) return p; } return NULL; } void insert(State* s){ states.push_back(s); } State* copy_state(State* s){ return new State(s); } State* find_minimized(State* s){ State* r = member(s); if(r == NULL){ r = copy_state(s); insert(r); } return r; } void set_transition(State* s, char c, State* t){ s->next[c] = t; } //辞書登録。inは辞書順ソートされている必要がある void create(const std::vector<std::string>& in, const std::vector<std::string>& out){ for(size_t i=0; i<MAX_WORD_SIZE; i++) TempState[i] = new State; PreviousWord = ""; TempState[0]->clear(); std::string CurrentWord = ""; std::string CurrentOutput = ""; for(size_t t=0; t<in.size(); t++){ CurrentWord = in[t]; CurrentOutput = out[t]; int i, j; i = 1; while(i<=CurrentWord.length() && i<=PreviousWord.length() && (CurrentWord[i-1] == PreviousWord[i-1])) i++; int PrefixLengthPlus1 = i; for(i=PreviousWord.size(); i>=PrefixLengthPlus1; i--){ set_transition(TempState[i-1], PreviousWord[i-1], find_minimized(TempState[i])); } for(i=PrefixLengthPlus1; i<=CurrentWord.length(); i++){ TempState[i]->clear(); set_transition(TempState[i-1], CurrentWord[i-1], TempState[i]); } if(in[t] != PreviousWord){ TempState[CurrentWord.length()]->final = true; } for(j=1; j<=PrefixLengthPlus1-1; j++){ if(TempState[j-1] == NULL) continue; std::string outputStr = TempState[j-1]->output[CurrentWord[j-1]]; std::string CommonPrefix = ""; for(int k=0; k<std::min(outputStr.length(), CurrentOutput.length()); k++){ if(outputStr[k] != CurrentOutput[k]) break; CommonPrefix += outputStr[k]; } std::string WordSuffix = outputStr.substr(CommonPrefix.length()); TempState[j-1]->output[CurrentWord[j-1]] = CommonPrefix; for(size_t c=0; c<0x100; c++){ if(TempState[j]->next[c] != NULL){ TempState[j]->output[c] = WordSuffix + TempState[j]->output[c]; } } CurrentOutput = CurrentOutput.substr(CommonPrefix.length()); } if(in[t] == PreviousWord){ TempState[CurrentWord.length()]->state_out += CurrentOutput; }else{ TempState[PrefixLengthPlus1-1]->output[CurrentWord[PrefixLengthPlus1-1]] = CurrentOutput; } PreviousWord = CurrentWord; } for(int i=CurrentWord.length(); i>=1; i--){ TempState[i-1]->next[PreviousWord[i-1]] = find_minimized(TempState[i]); InitialState = find_minimized(TempState[0]); } } //検索 std::string search(const std::string& query){ std::string ret = ""; State* p = InitialState; ret += p->state_out; for(int i=0; i<query.length(); i++){ if(!p->next[query[i]]) return ""; ret += p->output[query[i]]; p = p->next[query[i]]; ret += p->state_out; } if(p->final) return ret; return ""; } //内部状態を出力 void dump(){ std::queue<State*> que; que.push(InitialState); int state_id = 0; while(!que.empty()){ State* p = que.front(); que.pop(); if(p->_num >= 0) continue; p->_num = state_id; for(int i=0; i<0x100; i++){ if(p->next[i] == NULL) continue; que.push(p->next[i]); } state_id++; } std::cout << "size : " << states.size() << std::endl; std::map<int,int> memo; que.push(InitialState); while(!que.empty()){ State* p = que.front(); que.pop(); if(memo.count(p->_num) > 0) continue; memo[p->_num] = 1; std::cout << "node " << p->_num << " : [" << p->state_out << "] " << (p->final?"final":"") << std::endl; for(int i=0; i<0x100; i++){ if(p->next[i] == NULL) continue; std::cout << " " << (char)i << "->" << p->next[i]->_num << " [" << p->output[i] << "]" << std::endl; que.push(p->next[i]); } } } private: std::list<State*> states; State* TempState[MAX_WORD_SIZE]; std::string PreviousWord; State* InitialState; }; int main(){ Dictionary dic; std::vector<std::string> in, out; in.push_back("apr"); out.push_back("30"); in.push_back("aug"); out.push_back("31"); in.push_back("dec"); out.push_back("31"); in.push_back("feb"); out.push_back("28"); in.push_back("jan"); out.push_back("31"); in.push_back("jul"); out.push_back("31"); dic.create(in, out); dic.dump(); std::cout << "-------" << std::endl; std::cout << "apr : " << dic.search("apr") << std::endl; std::cout << "aug : " << dic.search("aug") << std::endl; std::cout << "dec : " << dic.search("dec") << std::endl; std::cout << "feb : " << dic.search("feb") << std::endl; std::cout << "jan : " << dic.search("jan") << std::endl; std::cout << "jul : " << dic.search("jul") << std::endl; std::cout << "abc : " << dic.search("abc") << std::endl; return 0; }
結果
size : 12 node 0 : [] a->1 [3] d->2 [31] f->3 [28] j->4 [31] node 1 : [] p->5 [0] u->6 [1] node 2 : [] e->7 [] node 3 : [] e->8 [] node 4 : [] a->9 [] u->10 [] node 5 : [] r->11 [] node 6 : [] g->11 [] node 7 : [] c->11 [] node 8 : [] b->11 [] node 9 : [] n->11 [] node 10 : [] l->11 [] node 11 : [] final ------- apr : 30 aug : 31 dec : 31 feb : 28 jan : 31 jul : 31 abc :
ランダムなアルファベット文字列でやってみても問題ないので一応大丈夫そう。
メモリ効率や高速化など実用レベルにするなら結構賢く書かないと大変そう。
GP-MIで遊ぶ
はじめに
http://live.nicovideo.jp/watch/lv228162988
先週のNL研のニコ生で、ベイズ的最適化についての招待講演を見ていて「SEは滑らかすぎる」という発言がよくわからなかったので、GP-MIを試してみる。
Contal et al., Gaussian Process Optimization with Mutual Information
http://arxiv.org/abs/1311.4825
http://econtal.perso.math.cnrs.fr/presentations/slides_icml14.pdf
コード
カーネルは、通常(？)のSquared Exponential kernel。
逆行列求めるところが重いのでEigen使う版も一緒に。
いつものごとく雑。
#include <iostream> #include <fstream> #include <sstream> #include <vector> #include <algorithm> #include <cmath> static const double PI = 3.14159265358979323846264338; //xorshift // 注意: longではなくint(32bit)にすべき unsigned long xor128(){ static unsigned long x=123456789, y=362436069, z=521288629, w=88675123; unsigned long t; t=(x^(x<<11)); x=y; y=z; z=w; return w=(w^(w>>19))^(t^(t>>8)); } //[0,1)の一様乱数 double frand(){ return xor128()%1000000/static_cast<double>(1000000); } //正規乱数 double normal_rand(double mu, double sigma2){ double sigma = sqrt(sigma2); double u1 = frand(), u2 = frand(); double z1 = sqrt(-2*log(u1)) * cos(2*PI*u2); //double z2 = sqrt(-2*log(u1)) * sin(2*PI*u2); return mu + sigma*z1; } #ifdef USE_EIGEN #include <Eigen/Core> #include <Eigen/LU> using namespace Eigen; #else //行列計算用 class Matrix { int m, n; std::vector<double> val; void splice(int i, int j){ int N = val.size(); if(i+(j-1)<N){ std::vector<double> vl(N-j, 0.0); int vl_idx = 0, val_idx = 0; while(val_idx<N){ if(val_idx<i||i+j-1<val_idx){ vl[vl_idx] = val[val_idx]; vl_idx++; } val_idx++; } val = vl; } } Matrix clone(){ Matrix ret(m, n); for(int i=0; i<m; i++){ for(int j=0; j<n; j++){ ret.setVal(i, j, val[n*i+j]); } } return ret; } public: Matrix():m(0),n(0){} Matrix(int m, int n):m(m),n(n){ if(m>0 && n>0){ for(int i=0; i<m*n; i++){ val.push_back(0.0); } } } Matrix(const Matrix& mat){ m = mat.m; n = mat.n; val = mat.val; } int getM() const { return m; } int getN() const { return n; } double getVal(int i, int j) const { return val[n*i+j]; } void setVal(int i, int j, double x){ val[n*i+j] = x; } bool isSquare(){ return n==m; } Matrix add(const Matrix& mat){ if(m == mat.m && n == mat.n){ Matrix ret(m, n); for(int i=0; i<m; i++){ for(int j=0; j<n; j++){ ret.setVal(i, j, val[n*i+j] + mat.getVal(i,j)); } } return ret; } return Matrix(-1, -1); } Matrix operator+(const Matrix& mat){ return add(mat); } Matrix sub(const Matrix& mat){ if(m == mat.m && n == mat.n){ Matrix ret(m, n); for(int i=0; i<m; i++){ for(int j=0; j<n; j++){ ret.setVal(i, j, val[n*i+j] - mat.getVal(i,j)); } } return ret; } return Matrix(-1, -1); } Matrix operator-(const Matrix& mat){ return sub(mat); } Matrix prod(const Matrix& mat){ if(n == 1 && mat.n == 1 && m == mat.m){ Matrix ret(m, 1); for(int i=0; i<m; i++){ ret.setVal(i, 0, getVal(i, 0) * mat.getVal(i, 0)); } return ret; } else if(n == mat.m){ Matrix ret(m, mat.n); for(int i=0; i<m; i++){ for(int j=0; j<mat.n; j++){ double d = 0.0; for(int k=0; k<n; k++){ d += val[n*i+k] * mat.getVal(k,j); } ret.setVal(i, j, d); } } return ret; } return Matrix(-1, -1); } Matrix operator*(const Matrix& mat){ return prod(mat); } void time(double x){ for(int i=0; i<m; i++){ for(int j=0; j<n; j++){ val[n*i+j] *= x; } } } Matrix transpose(){ Matrix ret(n, m); for(int i=0; i<m; i++){ for(int j=0; j<n; j++){ ret.setVal(j, i, val[n*i+j]); } } return ret; } double cofactor(int i, int j){ Matrix mat = clone(); mat.splice(i*mat.n, mat.m); mat.m -= 1; for(int k=mat.m-1; k>=0; k--){ mat.splice(k*mat.n+j, 1); } mat.n -= 1; return mat.det() * ( ((i+j+2)%2==0) ? 1 : -1); } double det(){ if(isSquare()){ if(m == 2){ return val[0]*val[3]-val[1]*val[2]; }else{ double d = 0; for(int k=0; k<n; k++){ d += val[k] * cofactor(0, k); } return d; } }else{ return 0.0; } } Matrix cofactorMatrix(){ Matrix mat(m, n); for(int i=0; i<m; i++){ for(int j=0; j<n; j++){ mat.setVal(j, i, cofactor(i, j)); } } return mat; } Matrix inverse(){ if(isSquare()){ double d = det(); if(d != 0){ Matrix mat; if(m>2){ mat = cofactorMatrix(); } else { mat = Matrix(2, 2); mat.setVal(0, 0, val[3]); mat.setVal(0, 1, -val[1]); mat.setVal(1, 0, -val[2]); mat.setVal(1, 1, val[0]); } mat.time(1 / d); return mat; } }else{ return Matrix(-1, -1); } return Matrix(-1, -1); } }; std::ostream& operator<<(std::ostream& os, const Matrix& mat){ for(int i=0; i<mat.getM(); i++){ for(int j=0; j<mat.getN(); j++){ os << mat.getVal(i,j) << " "; } os << std::endl; } return os; } #endif class GP_MI { double d; //xの次元 double alpha; //活用探索のトレードオフ調整用パラメータ double noise_sigma2; //ノイズの分散 int ITER; //サンプリング回数 double gamma; bool flg_obs; //推定,観測を繰り返しているかどうかの確認用フラグ //観測データ std::vector< std::vector<double> > obs_x; std::vector<double> obs_y; //各x_iの探索範囲 std::vector<double> x_min, x_max; //kernel(x,x')の値を求める double calc_kernel(const std::vector<double>& x1, const std::vector<double>& x2){ //Squared exponential kernel double sigma2 = 1.0; double l2 = 0.01; //比較的ガタガタを試すので小さめ double z = 0.0; for(size_t i=0; i<x1.size(); i++){ z += pow(((x1[i] - x2[i])), 2.0) / (2.0 * l2); } return sigma2 * exp(-z); } public: GP_MI(int d, double alpha, double noise_sigma2, int ITER): d(d), alpha(alpha), noise_sigma2(noise_sigma2), gamma(0.0), x_min(d, 0.0), x_max(d, 1.0), ITER(ITER), flg_obs(false) { } //xの探索範囲を指定 void set_x_range(int di, double mn, double mx){ x_min[di] = mn; x_max[di] = mx; } //現在までの観測結果から推定される関数をもとに、一番探した方がよいxを返す std::vector<double> get_best_x(){ //前回の観測が行われたかチェック if(!flg_obs){ std::cerr << "please observe and set" << std::endl; return std::vector<double>(); } flg_obs = false; //前計算C_T, Y_T #ifdef USE_EIGEN MatrixXd C_T = MatrixXd::Zero(obs_y.size(), obs_y.size()); MatrixXd Y_T = MatrixXd::Zero(obs_y.size(), 1); for(size_t i=0; i<obs_y.size(); i++){ for(size_t j=0; j<obs_y.size(); j++){ C_T(i,j) = calc_kernel(obs_x[i], obs_x[j]) + (i==j?noise_sigma2:0.0); } } for(size_t i=0; i<obs_y.size(); i++){ Y_T(i,0) = obs_y[i]; } MatrixXd C_T_inv = C_T.inverse(); #else Matrix C_T(obs_y.size(), obs_y.size()); Matrix Y_T(obs_y.size(), 1); for(size_t i=0; i<obs_y.size(); i++){ for(size_t j=0; j<obs_y.size(); j++){ C_T.setVal(i,j, calc_kernel(obs_x[i], obs_x[j]) + (i==j?noise_sigma2:0.0) ); } } for(size_t i=0; i<obs_y.size(); i++){ Y_T.setVal(i,0, obs_y[i]); } Matrix C_T_inv = C_T.inverse(); #endif //ITER回サンプリングしてargmaxを求める double best_a = -1.0e10; std::vector<double> best_x; double best_sigma2 = 0; for(size_t iter=0; iter<ITER; iter++){ //適当なxを生成 std::vector<double> x(d); for(size_t i=0; i<d; i++) x[i] = (x_max[i]-x_min[i])*frand() + x_min[i]; #ifdef USE_EIGEN MatrixXd kT = MatrixXd::Zero(obs_y.size(), 1); for(size_t i=0; i<obs_y.size(); i++){ kT(i,0) = calc_kernel(obs_x[i], x); } double kxx = calc_kernel(x, x); MatrixXd kCY = kT.transpose() * ( C_T_inv * Y_T ); MatrixXd kCk = kT.transpose() * ( C_T_inv * kT ); double mu = kCY(0,0); double sigma2 = kxx - kCk(0,0); #else Matrix kT(obs_y.size(), 1); for(size_t i=0; i<obs_y.size(); i++){ kT.setVal(i,0, calc_kernel(obs_x[i], x)); } double kxx = calc_kernel(x, x); Matrix kCY = kT.transpose() * ( C_T_inv * Y_T ); Matrix kCk = kT.transpose() * ( C_T_inv * kT ); double mu = kCY.getVal(0,0); double sigma2 = kxx - kCk.getVal(0,0); #endif //獲得関数の値を求める double phi = sqrt(alpha) * (sqrt(sigma2 + gamma) - sqrt(gamma)); if(best_a < mu + phi){ best_a = mu + phi; best_x = x; best_sigma2 = sigma2; } } gamma += best_sigma2; return best_x; } //get_best_x()で得られたxに対する実際の関数の値を、観測結果として登録 void set_observation(std::vector<double>& x, double y){ obs_x.push_back(x); obs_y.push_back(y); flg_obs = true; } }; //ブラックボックスな関数 double func(const std::vector<double>& x){ return x[0] * sin(10.0 * x[0]); } int main(){ //xの次元数 alpha noise_sigma2 サンプリング回数 GP_MI gpmi(1, 10.0, 0.0001, 500); //各次元の範囲 gpmi.set_x_range(0, 0.0, 1.0); //最初に何点か求めておく std::vector<double> x(1); x[0] = 0; gpmi.set_observation(x, func(x)); x[0] = 0.5; gpmi.set_observation(x, func(x)); x[0] = 1.0; gpmi.set_observation(x, func(x)); //自動で探索点を選ぶ for(int i=0; i<30; i++){ std::vector<double> best_x = gpmi.get_best_x(); std::cout << i << " : best_x = " << best_x[0] << std::endl; gpmi.set_observation(best_x, func(best_x)); } return 0; }
結果
最初に0,0.5,1.0の3点だけ求めた状態からスタート。
αがlog(2/δ), 0<δ<1とかだけど、小さいと隣の山に探しに行ってくれないので、大きいけど、α=10と50で試してみる。
グラフは、赤がμとσ^2、緑が獲得関数μ+φ。
x * sin(10.0 * x)
【α=10】
最適解は求められているけど、左の山が分散大きくても獲得関数の値が低いまま。
【α=50】
αを大きくする(大きくしすぎ？)と、左の山も探してくれてる。
x * cos(30.0 * x)
もっとガタガタしたもの。
【α=10】
うまくいっていない(本当の最適点はx=0.84ぐらい)。
動画で「滑らかすぎる」という発言があったけど、滑らかすぎて、関数の形の推定がうまくいかない＝適切な探索点をみつけられない、ということか(？)。
ただ、ARDなSEカーネルやMaternカーネルの実装、最適なパラメータ探索(MCMC)がちょっとよくわかってない・・・
http://www.slideshare.net/issei_sato/deep-learning-41310617
> BOツール公開予定
| ω・)待機。 |
参考
- http://live.nicovideo.jp/watch/lv228162988
- http://arxiv.org/abs/1311.4825
- http://econtal.perso.math.cnrs.fr/presentations/slides_icml14.pdf
- パターン認識と機械学習(下巻)
- Murphy, Machine Learning: a Probabilistic Perspective
Elman netを試す
はじめに
プロフェッショナルな「深層学習」本で紹介されているRNNの一種のElman netを書いてみる。
Recurrent Neural Network(RNN)とは
- 再帰型ニューラルネット
- ネットワーク内部に内部有向閉路を持つようなニューラルネットの総称
- Feedforwardの時は、入力層から出力層へ一方方向
- この構造のおかげで、時系列や言語モデル、系列ラベリングなど前の状態を持つような問題に対して考慮できる
- いろんな種類がある(以下はwikipediaから)
- Fully Recurrent network
- Hopfield network
- Boltzmann machine
- Simple recurrent network
- Elman net
- Jordan net
- Echo state network
- Long short term memory(LSTM) network
- Bi-directional RNN
- Hierarchical RNN
- Stochastic neural networks
- Fully Recurrent network
Simple recurrent network
- 言葉のとおり、シンプルな構造のRNN
- 代表的には「Elman net」と「Jordan net」がある
- 3層(入力層、隠れ層、出力層)のFeedforward Neural netについて、どこからどこへ閉路を作るかで異なる
- Elman net : 隠れ層から隠れ層への辺がある
- Jordan net : 出力層から隠れ層への辺がある
コード
「深層学習」本では、Elman netのBPTT法(Backpropagation through time、時間でネットワークを分けて普通に逆誤差伝播)での学習が紹介されていると思われるので、これを書いてみる。
#include <iostream> #include <vector> #include <cstdio> #include <algorithm> #include <cmath> static const double PI = 3.14159265358979323846264338; //xorshift // 注意: longではなくint(32bit)にすべき unsigned long xor128(){ static unsigned long x=123456789, y=362436069, z=521288629, w=88675123; unsigned long t; t=(x^(x<<11)); x=y; y=z; z=w; return w=(w^(w>>19))^(t^(t>>8)); } //[0,1)の一様乱数 double frand(){ return xor128()%1000000/static_cast<double>(1000000); } //正規乱数 double normal_rand(double mu, double sigma2){ double sigma = sqrt(sigma2); double u1 = frand(), u2 = frand(); double z1 = sqrt(-2*log(u1)) * cos(2*PI*u2); //double z2 = sqrt(-2*log(u1)) * sin(2*PI*u2); return mu + sigma*z1; } //隠れ層が1層の単純なRNN // BPTT(Backpropagation through time)法で学習 class SimpleRNN { double eps; //学習率 int Nin, Nhid, Nout; //各層のユニット数(バイアスを除く) std::vector< std::vector<double> > W; //時刻t-1での隠れ層から時刻tでの隠れ層へ std::vector< std::vector<double> > Win; //時刻tでの入力層から隠れ層へ std::vector< std::vector<double> > Wout; //時刻tでの隠れ層から出力層へ std::vector< std::vector<double> > u, v; //各時刻での隠れ層の入力, 出力層の入力 std::vector< std::vector<double> > z, y; //各時刻での隠れ層の出力, 出力層の出力 std::vector< std::vector<double> > delta, delta_out; //各時刻での各ユニットのδ値 public: SimpleRNN(int Nin, int Nhid, int Nout, double eps): eps(eps), Nin(Nin), Nhid(Nhid), Nout(Nout), W(Nhid+1, std::vector<double>(Nhid, 0)), Win(Nin+1, std::vector<double>(Nhid, 0)), Wout(Nhid+1, std::vector<double>(Nout, 0)) { for(int i=0; i<Nhid+1; i++){ for(int j=0; j<Nhid; j++){ W[i][j] = normal_rand(0.0, 0.1); } } for(int i=0; i<Nin+1; i++){ for(int j=0; j<Nhid; j++){ Win[i][j] = normal_rand(0.0, 0.1); } } for(int i=0; i<Nhid+1; i++){ for(int j=0; j<Nout; j++){ Wout[i][j] = normal_rand(0.0, 0.1); } } } std::vector<int> forward_propagation(const std::vector< std::vector<double> >& in){ int T = in.size(); u = std::vector< std::vector<double> >(T, std::vector<double>(Nhid, 0)); v = std::vector< std::vector<double> >(T, std::vector<double>(Nout, 0)); z = std::vector< std::vector<double> >(T, std::vector<double>(Nhid+1, 0)); y = std::vector< std::vector<double> >(T, std::vector<double>(Nout, 0)); //各時刻tでの値を求める std::vector<int> ret; for(int t=0; t<T; t++){ for(int i=0; i<Nhid; i++){ u[t][i] = 0.0; } //入力層の出力->隠れ層の入力 for(int i=0; i<Nin; i++){ for(int j=0; j<Nhid; j++){ u[t][j] += in[t][i] * Win[i][j]; } } for(int i=0; i<Nhid; i++){ //バイアス u[t][i] += 1.0 * Win[Nin][i]; } //時刻t-1での隠れ層の出力->時刻tでの隠れ層の入力 for(int i=0; i<Nhid+1; i++){ if(t!=0){ for(int j=0; j<Nhid; j++){ u[t][j] += z[t-1][i] * W[i][j]; } } if(t==0 && i==Nhid){ //バイアス for(int j=0; j<Nhid; j++){ u[t][j] += 1.0 * W[Nin][j]; } } } //時刻tでの隠れ層の出力 for(int i=0; i<Nhid; i++){ z[t][i] = 1.0 / (1.0 + exp(-u[t][i])); } z[t][Nhid] = 1.0; //時刻tでの出力層の入力 for(int i=0; i<Nhid+1; i++){ for(int j=0; j<Nout; j++){ v[t][j] += z[t][i] * Wout[i][j]; } } //時刻tでの出力層の出力 double Z = 0; double mx = -1.0; int mx_i = -1; for(int i=0; i<Nout; i++){ Z += exp(v[t][i]); } for(int i=0; i<Nout; i++){ y[t][i] = exp(v[t][i]) / Z; if(mx < y[t][i]){ mx = y[t][i]; mx_i = i; } } ret.push_back(mx_i); } return ret; } double back_propagation(const std::vector< std::vector<double> >& in, const std::vector<int>& out){ double err = 0.0; int T = in.size(); std::vector<int> res = forward_propagation(in); delta = std::vector< std::vector<double> >(T, std::vector<double>(Nhid, 0)); delta_out = std::vector< std::vector<double> >(T, std::vector<double>(Nout, 0)); //時刻の逆方向からdelta値と重みの更新を連続的にしていく for(int t=T-1; t>=0; t--){ //時刻tでの出力層のdelta値 for(int i=0; i<Nout; i++){ if(out[t] == i){ delta_out[t][i] = y[t][i] - 1.0; err += - 1.0 * log(y[t][i]); }else{ delta_out[t][i] = y[t][i] - 0.0; } } //時刻tでの隠れ層の出力層からのdelta値 for(int i=0; i<Nhid; i++){ for(int j=0; j<Nout; j++){ double fu = 1.0 / (1.0 + exp(-u[t][i])); double fdash = fu * (1.0 - fu); delta[t][i] += Wout[i][j] * delta_out[t][j] * fdash; } } //時刻tでの隠れ層の時刻t+1での隠れ層からのdelta値 for(int i=0; i<Nhid; i++){ for(int j=0; j<Nhid; j++){ if(t+1>=T) continue; //時刻T以降の場合はdelta=0として扱う double fu = 1.0 / (1.0 + exp(-u[t][i])); double fdash = fu * (1.0 - fu); delta[t][i] += W[i][j] * delta[t+1][j] * fdash; } } //delta値を使って各重みの微分を計算し、SGDで重みを更新 for(int j=0; j<Nhid+1; j++){ //隠れ層->出力層 for(int k=0; k<Nout; k++){ Wout[j][k] -= eps * (delta_out[t][k] * z[t][j]); } } for(int pj=0; pj<Nhid+1; pj++){ //隠れ層->隠れ層 for(int j=0; j<Nhid; j++){ if(t-1>=0){ W[pj][j] -= eps * (delta[t][j] * z[t-1][pj]); }else{ //時刻t<0の場合は、バイアスだけ1でそれ以外の出力は0として扱う if(pj==Nhid){ //バイアス W[pj][j] -= eps * (delta[t][j] * 1.0); } } } } for(int i=0; i<Nin+1; i++){ //入力層->隠れ層 for(int j=0; j<Nhid; j++){ if(i<Nin){ Win[i][j] -= eps * (delta[t][j] * in[t][i]); }else{ //バイアス Win[i][j] -= eps * (delta[t][j] * 1.0); } } } } return err; } }; int main(){ const int TERM = 500; //TERM個単位でerrを確認 const double finish_err = 0.01; //errがfinish_err未満になったら学習終了 int Nin, Nhid, Nout; double eps; int trainN, testN; std::cin >> Nin >> Nhid >> Nout; std::cin >> eps; /// 学習 //////////////////// std::cin >> trainN; std::vector< std::vector< std::vector<double> > > in; std::vector< std::vector<int> > out; for(int i=0; i<trainN; i++){ int t; std::cin >> t; std::vector< std::vector<double> > in_one; std::vector<int> out_one; for(int j=0; j<t; j++){ double val; int outval; std::cin >> outval; std::vector<double> v; for(int k=0; k<Nin; k++){ std::cin >> val; v.push_back(val); } in_one.push_back(v); out_one.push_back(outval); } in.push_back(in_one); out.push_back(out_one); } //データシャッフル用 std::vector<int> rnd; for(int i=0; i<trainN; i++){ rnd.push_back(i); } std::random_shuffle(rnd.begin(), rnd.end()); //学習ループ SimpleRNN rnn(Nin, Nhid, Nout, eps); int iter = 0; double err = 0.0; while(true){ //TERM個単位で確認 if(iter!=0 && iter%TERM == 0){ err /= TERM; std::cerr << "err = " << err << std::endl; if(err < finish_err) break; err = 0; } err += rnn.back_propagation(in[rnd[iter%trainN]], out[rnd[iter%trainN]]); iter++; } /// テスト //////////////////// std::cin >> testN; for(int i=0; i<testN; i++){ int t; std::cin >> t; std::vector< std::vector<double> > in_one; std::vector<int> out_one; for(int j=0; j<t; j++){ double val; int outval; std::cin >> outval; std::vector<double> v; for(int k=0; k<Nin; k++){ std::cin >> val; v.push_back(val); } in_one.push_back(v); out_one.push_back(outval); } //結果の出力 std::vector<int> res = rnn.forward_propagation(in_one); std::cout << "ref:"; for(int j=0; j<out_one.size(); j++){ std::cout << " " << out_one[j]; } std::cout << std::endl; std::cout << "out:"; for(int j=0; j<res.size(); j++){ std::cout << " " << res[j]; } std::cout << std::endl; std::cout << std::endl; } return 0; }
入力ファイルの形式
[入力層のユニット数] [隠れ層のユニット数] [出力層のユニット数] [学習率] [学習インスタンス数] [ケース1の数] [出力1] [次元1の値] [次元2の値] ... [出力2] [次元1の値] [次元2の値] ... [出力3] [次元1の値] [次元2の値] ... ... [テストインスタンス数] [ケース1の数] [出力1] [次元1の値] [次元2の値] ... [出力2] [次元1の値] [次元2の値] ... [出力3] [次元1の値] [次元2の値] ... ...
結果
おもちゃなケースで確認だけ。
http://d.hatena.ne.jp/jetbead/20111125/1322236626 でのサンプルを今回の形式に変換したもの。
- 入力
10 10 3 0.01 6 3 2 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 2 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 2 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 4 2 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 2 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 1 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 3 2 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 2 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 2 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 4 2 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 2 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 1 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 3 2 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 2 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 2 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 4 2 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 2 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 2 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 2 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 6 3 2 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 2 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 2 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 4 2 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 2 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 1 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 3 2 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 2 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 2 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 4 2 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 2 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 1 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 3 2 1.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 2 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 0.0 2 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 4 2 0.0 0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 2 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 0.0 2 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0 0.0 2 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 1.0
- 出力
$ ./a.out < in err = 2.20794 err = 1.84924 err = 1.54379 err = 1.22335 err = 0.931081 err = 0.681177 err = 0.493641 err = 0.360733 err = 0.265542 err = 0.198631 err = 0.153081 err = 0.121373 err = 0.0987013 err = 0.0824887 err = 0.0702189 err = 0.0606353 err = 0.0533114 err = 0.0473444 err = 0.0423609 err = 0.0383958 err = 0.034996 err = 0.0320187 err = 0.0295968 err = 0.0274426 err = 0.0254882 err = 0.0238809 err = 0.0224107 err = 0.0210389 err = 0.0199064 err = 0.0188471 err = 0.0178354 err = 0.0170008 err = 0.0162053 err = 0.0154302 err = 0.0147935 err = 0.0141766 err = 0.0135647 err = 0.0130652 err = 0.0125743 err = 0.0120795 err = 0.0116788 err = 0.0112797 err = 0.0108714 err = 0.010544 err = 0.0102138 err = 0.00987132 ref: 2 2 2 out: 2 2 2 ref: 2 2 0 1 out: 2 2 0 1 ref: 2 2 2 out: 2 2 2 ref: 2 2 0 1 out: 2 2 0 1 ref: 2 2 2 out: 2 2 2 ref: 2 2 2 2 out: 2 2 2 2
ちゃんと収束して学習データはちゃんとラベルつけられている。
しかし、展開する時間tが長くなるほどネットワークが深く(deep)なってしまい、勾配消失(または発散)問題が起きてしまってうまく学習できない。
参考
- 岡谷, 機械学習プロフェッショナルシリーズ「深層学習」, 講談社
- http://www.slideshare.net/beam2d/pfi-seminar-20141030rnn
- http://qiita.com/icoxfog417/items/2791ee878deee0d0fd9c
- https://en.wikipedia.org/wiki/Types_of_artificial_neural_networks#Recurrent_neural_network
- https://github.com/mattya/RNN-colle/wiki/Tutorial_jp
- http://wbawakate.jp/wp-content/uploads/2015/03/RNN%E3%83%95%E3%82%9A%E3%83%AC%E3%82%BB%E3%82%99%E3%83%B3.pdf
多層ニューラルネットを試す
はじめに
FeedForwardNeuralNetwork。プロフェッショナルな「深層学習」本のバックプロパゲーションの導出が丁寧にされていてわかりやすかったので、それに合わせて書いてみる。
各層の活性化関数はロジスティック(シグモイド)関数、出力層の活性化関数はソフトマックス関数、誤差関数は交差エントロピー。
コード
1インスタンスごとに重みを更新(SGD)。
直近TERM個のインスタンスの誤差平均が終了条件を満たしたら終了。
学習が収束してくれているので大丈夫そう。
#include <iostream> #include <vector> #include <cstdio> #include <algorithm> #include <cmath> static const double PI = 3.14159265358979323846264338; //xorshift // 注意: longではなくint(32bit)にすべき unsigned long xor128(){ static unsigned long x=123456789, y=362436069, z=521288629, w=88675123; unsigned long t; t=(x^(x<<11)); x=y; y=z; z=w; return w=(w^(w>>19))^(t^(t>>8)); } //[0,1)の一様乱数 double frand(){ return xor128()%1000000/static_cast<double>(1000000); } //正規乱数 double normal_rand(double mu, double sigma2){ double sigma = sqrt(sigma2); double u1 = frand(), u2 = frand(); double z1 = sqrt(-2*log(u1)) * cos(2*PI*u2); //double z2 = sqrt(-2*log(u1)) * sin(2*PI*u2); return mu + sigma*z1; } struct UnitW { double u; //ユニットへの入力 double z; //ユニットの出力 double delta; //ユニットでのδ_j^l = ∂E_n / ∂u_j^(l) std::vector<double> w; //このユニットから次の層のユニットへの重み UnitW(int num_next_layer_node):w(num_next_layer_node, 0){ for(int i=0; i<w.size(); i++){ w[i] = normal_rand(0.0, 0.1); } } }; class MultiLayerNeuralNetwork { double eps; //学習率 std::vector< std::vector<UnitW> > network; void init_network(const std::vector<int>& num_layer_node){ for(int i=0; i<num_layer_node.size(); i++){ network.push_back( std::vector<UnitW>() ); //(num_layer_node[i])番目をバイアスとして使うのですべての層で+1用意 for(int j=0; j<num_layer_node[i] + 1; j++){ if(i < num_layer_node.size()-1){ //入力層、隠れ層 network[i].push_back(UnitW(num_layer_node[i+1])); if(j == num_layer_node[i]){ //バイアスの重みは0 for(int k=0; k<network[i][j].w.size(); k++){ network[i][j].w[k] = 0.0; } } }else{ //出力層 network[num_layer_node.size()-1].push_back(UnitW(0)); } } } } public: MultiLayerNeuralNetwork(double eps, const std::vector<int>& num_layer_node):eps(eps){ init_network(num_layer_node); } std::vector<double> forward_propagation(const std::vector<double>& in){ for(int i=0; i<network.size(); i++){ for(int j=0; j<network[i].size(); j++){ network[i][j].u = 0.0; } } //入力層 for(int j=0; j<network[0].size()-1; j++){ network[0][j].z = in[j]; } network[0][network[0].size()-1].z = 1.0; //隠れ層と出力層 for(int i=0; i<network.size()-1; i++){ for(int j=0; j<network[i].size(); j++){ //現在の層の出力zと重みを合わせて次の層の入力に渡す for(int k=0; k<network[i+1].size()-1; k++){ network[i+1][k].u += network[i][j].z * network[i][j].w[k]; } } //次の層の出力zを計算 if(i < network.size()-2){ //隠れ層 for(int j=0; j<network[i+1].size(); j++){ network[i+1][j].z = 1.0 / (1.0 + exp(-network[i+1][j].u)); } network[i+1][network[i+1].size()-1].z = 1.0; //バイアスは常に出力は1 }else{ //出力層 double Z = 0.0; for(int j=0; j<network[i+1].size()-1; j++){ Z += exp(network[i+1][j].u); } for(int j=0; j<network[i+1].size()-1; j++){ network[i+1][j].z = exp(network[i+1][j].u) / Z; } } } std::vector<double> ret; for(int i=0; i<network[network.size()-1].size()-1; i++){ ret.push_back(network[network.size()-1][i].z); } return ret; } double back_propagation(const std::vector<double>& in, const std::vector<double>& d){ double err = 0.0; for(int i=0; i<network.size(); i++){ for(int j=0; j<network[i].size(); j++){ network[i][j].delta = 0.0; } } //現在のネットワークでの結果を計算 std::vector<double> y = forward_propagation(in); //出力層のdelta値の計算 for(int i=0; i<network[network.size()-1].size()-1; i++){ network[network.size()-1][i].delta = y[i] - d[i]; err += - d[i] * log(y[i]); } //各ユニットのdelta値の計算 for(int i=network.size()-2; i>=0; i--){ for(int j=0; j<network[i].size(); j++){ for(int k=0; k<network[i+1].size()-1; k++){ double fu = 1.0 / (1.0 + exp(-network[i][j].u)); double fdash = fu * (1.0 - fu); network[i][j].delta += network[i+1][k].delta * network[i][j].w[k] * fdash; } } } //各ユニットのdelta値を使って各重みにおける微分を計算し、SGDで重みを更新 for(int i=network.size()-2; i>=0; i--){ for(int j=0; j<network[i].size(); j++){ for(int k=0; k<network[i+1].size()-1; k++){ network[i][j].w[k] -= eps * (network[i+1][k].delta * network[i][j].z); } } } return err; } }; int main(){ const int TERM = 500; //TERM個単位でerrを確認 const double finish_err = 0.01; //errがfinish_err未満になったら学習終了 int layerN, trainN, testN; double eps; std::cin >> layerN; //各層のユニット数 std::vector<int> num_layer_node; for(int i=0; i<layerN; i++){ int num; std::cin >> num; num_layer_node.push_back(num); } int inN = num_layer_node[0]; int outN = num_layer_node[num_layer_node.size()-1]; std::cin >> eps; /// 学習 //////////////////// std::cin >> trainN; std::vector< std::vector<double> > out, in; for(int i=0; i<trainN; i++){ std::cerr << "\r" << i; std::vector<double> out_one(outN, 0), in_one; int t; std::cin >> t; out_one[t] = 1; out.push_back(out_one); for(int j=0; j<inN; j++){ double d; std::cin >> d; in_one.push_back(d); } in.push_back(in_one); } std::cerr << std::endl; //データシャッフル用 std::vector<int> rnd; for(int i=0; i<trainN; i++){ rnd.push_back(i); } std::random_shuffle(rnd.begin(), rnd.end()); //学習ループ MultiLayerNeuralNetwork nn(eps, num_layer_node); int iter = 0; double err = 0.0; while(true){ //TERM個単位で確認 if(iter != 0 && iter%TERM == 0){ err /= TERM; std::cerr << "err = " << err << std::endl; if(err < finish_err) break; err = 0; } err += nn.back_propagation(in[rnd[iter%trainN]], out[rnd[iter%trainN]]); iter++; } /// テスト //////////////////// std::cin >> testN; int match = 0; for(int i=0; i<testN; i++){ std::vector<double> in_one; int out_t; std::cin >> out_t; for(int j=0; j<inN; j++){ double d; std::cin >> d; in_one.push_back(d); } std::vector<double> res = nn.forward_propagation(in_one); double res_mx = 0; int res_i = 0; for(int j=0; j<res.size(); j++){ if(res[j]>res_mx){ res_mx = res[j]; res_i = j; } } if(out_t == res_i) match++; std::cout << out_t << "\t" << res_i << "\t" << res_mx << std::endl; } //Accuracyの出力 std::cout << "Acc = " << match/(double)testN << std::endl; return 0; }
入力ファイルの形式
[layer数] [入力層のユニット数] [隠れ層1層目のユニット数] ... [隠れ層k層目のユニット数] [出力層のユニット数] [学習率] [学習インスタンス数] [出力] [次元1の値] [次元2の値] ... ... テストインスタンス数 [出力] [次元1の値] [次元2の値] ... ...
出力は、0〜出力ユニット数-1の値。入力の次元は、入力層のユニット数と同じ事が必要。
結果
隠れ層1層のケースで確認。
XORパターン
- 入力
$ cat xor.txt 3 2 2 2 0.01 4 0 0 0 1 0 1 1 1 0 0 1 1 4 0 0 0 1 0 1 1 1 0 0 1 1
- 出力
$ ./a.out < xor.txt 3 err = 0.697136 ... err = 0.00992856 0 0 0.992468 1 1 0.991095 1 1 0.991083 0 0 0.986047 Acc = 1
MNISTの手書き数字認識
scikit_learn_dataのやつ。
- 入力
$ cat mnist.txt 3 784 100 10 0.01 60000 0 0.0 0.0 ...(省略) ...(省略) 10000 0 0.0 0.0 ...(省略) ...(省略)
- 出力
$ ./a.out < mnist.txt 59999 err = 2.28309 err = 1.34394 err = 1.08896 ... Acc = 0.9776
- 収束までにかかる時間
CPU(i7-4790)計算で、結果でるまで16分30秒ぐらい。
コードが違うけど、chainerだとGPU使って30秒程度で終わる、、、やっぱりchainerだな！
参考
- 岡谷, 機械学習プロフェッショナルシリーズ「深層学習」, 講談社
- http://aidiary.hatenablog.com/entry/20140205/1391601418
- 以前遊んだやつは、http://d.hatena.ne.jp/jetbead/20121110/1352525897
Feature-Weighted Linear Stackingメモ
はじめに
Sill et al., Feature-Weighted Linear Stacking, 2009
http://arxiv.org/abs/0911.0460
最近、コンペ上位者の手法としてよく見かける手法「Stacking」の一つについてメモ。
Stacking
Feature-Weighted Linear Stacking(FWLS)
概要的には、各学習結果をまとめあげる部分を線形和として表す際に、メタ素性をさらに考慮することで、細かい調整を行えるようにした感じ。
- L個の機械学習モデル g_i
- 入力xに対して、g_i(x)は予測結果を表す
- M個のメタ素性関数 f_i
- 入力xに対して、f_i(x)はメタ素性に関する何らかの実数を表す
Standard linear regression stacking
Feature-weighted linear stacking
FWLS optimization problem
- 最適化したい問題(v_{ij}の学習)は、
- min_v Σ_入力x Σ_{i,j} (v_{ij} f_j(x) g_i(x) - 正解y(x))^2
- ただし、入力xはスタッキングのパラメータを求めるために各学習器の学習で使っていない入力集合を使ったほうがよい
- K-fold cross validationとか、各学習に使う入力集合には気をつける
メタ素性の例
参考
- http://arxiv.org/abs/0911.0460
- http://ibisforest.org/index.php?%E3%83%A1%E3%82%BF%E5%AD%A6%E7%BF%92
- http://shindannin.hatenadiary.com/entry/2015/03/13/101945
- http://www.chioka.in/stacking-blending-and-stacked-generalization/
- http://en.wikipedia.org/wiki/Ensemble_learning#Stacking
- http://en.wikipedia.org/wiki/Meta_learning_(computer_science)
- 杉山ら, 統計的学習の基礎 データマイニング・推論・予測, 共立出版
- Murphy, Machine Learning A Probablistic Perspective
- Aggarwal, Data Classification Algorithms and Applications
- (下2つの本ではこの論文はReferenceされている)
係り受け解析メモ
はじめに
今年の目標にしていた係り受け解析関係の資料について雑多にメモしておく。リンク集。
拾いきれていない、最新の論文まで追えていないので、あとで追加・整理し直す。
文節単位がよいか、単語単位がよいかの議論
解析処理について説明している日本語資料
- 海野, 統計的係り受け解析入門
- G. Neubig, 自然言語処理プログラミング勉強会12 係り受け解析
- 黒橋, 結構やるな, KNP
- 笹野, 構文・述語項構造解析システムKNPの解析の流れと特徴
- 海野さんの論文紹介
- http://blog.unnono.net/2010/01/nlp20102-non-projective-dependency.html
- McDonaldら, Non-Projective Dependency Parsing using Spanning Tree Algorithm, HLT-EMNLP2005
- Smithら, Probabilistic Models of Nonprojective Dependency Trees, EMNLP2007
- Martinsら, Concise Integer Linear Programming Formulations for Dependency Parsing, ACL2009
- 颯々野, 日本語係り受け解析の線形時間アルゴリズム
- 工藤ら, Support Vector Machineによる日本語係り受け解析
- 春野ら, 決定木を用いた日本語係り受け解析
- 能地, 最近のTransition-based Parsing, PFI seminar
その他の話題
- Kooら, Simple Semi-supervised Dependency Parsing
- 教師なし
- 双対分解の係り受け解析への適用
- 海野ら, N-gram統計量からの係り受け情報の復元
解析器
- KNP
- CaboCha
- J.DepP
- EDA(Easily adaptable Dependency Analyzer)係り受け解析器
- MaltParser
- MSTParser(Minimum-Spanning Tree Parser)
- Stanford Parser
係り受け解析について説明している書籍
ネット上の文書をきれいにするライブラリ「NetTextCleaner」を公開しました
エイプリルフールネタです:) 見ればわかりますが、30分ぐらいで作ったものなので、言ってることもコードもデータも結構適当です:p
はじめに
インターネット上で扱われる文書(Twitterや2ch、ニコニコ動画など)には、特殊な用語や言い回しをしているものがあり、ネットの文書を扱う多くの自然言語処理研究者を苦しめています。
一番問題となることは、「これらの特殊な用語や言い回しを含む文を形態素解析すると解析に失敗すること」で、これによってIDFの計算や統計的言語モデルがおかしくなったりしてしまいます。
そこで、「書き手が書きたいことが、ひねくれて、実際に書かれる」という仮説を考え、生成的な確率モデルである「Noisy Channel Model」を用いてモデル化し、書き手が本当に書きたかったきれいな文へ変換することで、自然言語処理に失敗しない文章にすることができるライブラリを作成しました。
今後、各ウェブ企業などに導入され、ネット上のあらゆる文書に適用されることで、きれいになっていくでしょう。
Noisy Channel Modelによる変換
本当に書きたかったきれいな文章をS、実際に書いてしまった特殊な言い回しをした文章をTとします。
これをNoisy Channel Modelを用いると、
のようにモデル化でき、変換モデルを統計的に求めることができる気がします。
本ライブラリでは、ネット上からサンプルを各1件ずつ抽出し、統計的に変換モデルを作成しています。
ライブラリの概要
「変換モデル」と「変換プログラム」の2つから成ります。
「変換プログラム」は、正規表現を用いて上記の変換システムを実装してあり、「変換モデル」は変換モデルが含まれています。
実際のプログラムは、下記のAppendixにまとめられています。
使用例
以下のように、ネット上の文書を、一般的な理解可能な文書へ変換します。
$ ls model.dat cleaner.pl $ perl cleaner.pl あくしろよ => はやくしてください 猫、ぐうかわ => 猫、ぐうの音もでないほどかわいい 賢者タイム => 男性のみに訪れる特別な時間 帰って、どうぞ => 帰ってください
評価
Twitterのタイムラインを流れていたツイートに対してこれらを適用した結果、精度66.7%(3件中2件正解)という高精度で変換することができています。
免責事項
このライブラリを用いて起きたあらゆる問題・障害に対して責任を持ちません。
おそらく適当に変更して使っていただいて問題ないと思います。
Appendix
変換プログラムと変換モデルを以下に添付します。コピペによって利用する事ができます。
【cleaner.pl】
#!/usr/local/bin/perl use strict; use warnings; my $dic_path = "./model.dat"; my @src_rules; my @dst_rules; open(IN, "<", $dic_path) or die; while(<IN>){ chomp; my ($src_rule, $dst_rule) = split(/\s+/, $_, 2); push(@src_rules, $src_rule); push(@dst_rules, $dst_rule); } while(<>){ chomp; my $text = $_; for(my $i=0; $i<@src_rules; $i++){ $text =~ s/$src_rules[$i]/qq(\"$dst_rules[$i]\")/eeg; } print "=> ", $text, "\n"; }
【model.dat】
(.+?)かわいいよ\1 $1がかわいい (.+?)\1アンド\1 とても$1 (.+?)のバーゲンセール $1の希少価値が落ちている (.+?)ェ… $1… (.+?)過ぎて謝謝 $1過ぎて申し訳ない (.+?)に(\d+?)ペリカ $1だと思います (.+?)に自信ニキ $1に自信がある人 (.+?)は(ポイー|ﾎﾟｲ-)で $1はいらない あ…(察し) 心中お察しします アイエ.+?\!\s*(.+?)\!\?\s*\1ナンデ\!\? $1だ あくしろよ はやくしてください アッー！ アーッ！ あまり私を怒らせない方がいい 私を怒らせると大変な事になりますよ ありがとナス ありがとう ＼(.+?)だー！／ $1です あーげない あげません いーらない いりません くーださい ください ｲｲﾊﾅｼﾀﾞﾅｰ いい話です 1ヶ月綺麗な言葉をかけ続けた(.+?)と1ヶ月罵声を浴びせ続けた\1 1ヶ月きれいな言葉と罵声で$1がどう変化したか いつから(.+?)だと錯覚していた？ $1は錯覚です 今北産業 今来たところなのでこれまでの流れを三行で説明してください びゃあ゛ぁ゛゛ぁうまひぃ゛ぃぃ゛ いやあ、うまい (.+?)しろよ、色々と捗るぞ $1すれば世界が変わりますよ？オススメですよ？ インド人を右に ハンドルを右に イﾞェアアアア \(断末魔\) (.+?)から来るぞ！気をつけろ！ $1から何か来ます ウキウキでワロタｗ+? ウキウキで笑いました イライラでワロタｗ+? イライラで笑いました (?<!ドワ|ニワ)ンゴ 。 (.+?)の法則が乱れる！ (.+?)がよくわからないがとにかくめちゃくちゃな状態になっている うはｗｗｗおｋｗｗｗ ハイテンションなので笑って許せます うーんこの なにかがおかしい 教えてエロい人 教えてください オナシャス お願いします おめぇの(.+?)ねぇから あなたの$1はありません おもんねーわ 面白くない オワコン 終わったコンテンツ (.+?)がね・・・ $1が致命的です (.+?)か何か？ $1ではないと思いますが$1でしょうか？ キマシタワー 百合展開、きました きんもーっ☆ キモいです は犠牲になったのだ は死にました ぎゃんかわ とてもかわいい ギャン泣き ギャンギャン泣きわめく 草不可避 笑いを耐えられません くやしいのう くやしいですね ぐうかわ ぐうの音もでないほどかわいい グンマー 群馬県 賢者タイム 男性のみに訪れる特別な時間 激おこぷんぷん丸 結構怒っています こ↑こ↓ ここ これで勝つる これで勝てます GOD’S and DEATH ありがとうございます 33\-4 圧倒的大差 (.+?)、どうぞ $1ください 志村うしろ！ 後ろから危険が迫っています 正直、すまんかった ごめんなさい （小並感） （小学生並みの感想） が一晩でやってくれました。 が尋常ではないはやさで仕事をしてくれました じゃあの また会いましょう 少し、頭冷やそうか 冷静になりましょう すっごい(.+?)よ！ すごく$1です！ すまぬ…すまぬ… ごめんなさい ですしおすし ですし でっていう だからなに？ どう見ても(.+?)です。本当にありがとうございました。 どう見ても$1です\(残念\) 日本語でおk 日本語でお願いします 日本よ、これが(.+?)だ これこそ真の$1 はいじゃないが それどころじゃないでしょう 微レ存 微粒子レベルで存在しています ふつくしい 美しい △ さんかっこいい ボスケテ 助けて 無理ぽ もう無理っぽい メシウマ 他人の不幸で今日も飯がうまい メシマズ 他人の幸福で飯がまずい もうだめぽ もうダメだ・・・ (.+?)のライフはゼロよ！ $1は瀕死です！ モルダーあなた疲れてるのよ 疲れておかしくなっています やめたげてよぉ やめてあげてください やろう、ぶっころしてやる ブチ切れました 許してヒヤシンス 許してください よろしい、ならば戦争だ でしたら、戦争で勝負をつけましょう ワシの(.+?)は108式まであるぞ まだ本気を出していないが、本気を出したらもっとすごいですよ？
