2014-08-24

Reservoir sampling

$N$ 件のデータから $k$ 件をランダムサンプリングする場合のアルゴリズムで、
$N$ が未知数であっても、 $N$ を求める必要がなく、
1passでランダムサンプリングができる。
今まで $N$ を求めるためにかかっていた時間を家事や育児にあてることができる。

概要

$N$ 件のデータを $S[0],S[1], \cdots ,S[N-1]$ とする。
このアルゴリズムでは、Reservoir（貯水池）を $R[0],R[1], \cdots ,R[k-1]$ として、
そこに $k$ 件のデータを随時貯めておく。
Reservoirは処理の過程でどんどん更新されていき、
最終的にはこれがサンプリング結果となる。

アルゴリズム

アルゴリズムは基本的に２つのステップから構成される。

ステップ①
最初の $k$ 件をとりあえずReservoirにぶっこんでおく。
$\displaystyle R[i ] \left( 0\leq i < k \right)$ にデータ $S[i]$ を代入。

ステップ②
データ $S[i] \left( k \leq i < N \right)$ については $\frac{k}{i+1}$ の確率でReservoirの保持内容を更新する。
$[0,i]$ の整数値をとる一様乱数 $j$ を生成し、 $j < k$ の場合、 $R[j]$ に $S[i]$ を代入。

解釈

本当にランダムサンプリングされているのか気になるところ。
２つの場合に分けてデータ $S [ i ]$ が選ばれる確率を考えてみる。

データ $S [ i ] \left( i < k \right)$ がReservoirに残る確率
$\hspace{3em} = R [ i ]$ が更新されない確率
$\hspace{3em} = \frac{k}{k+1} \frac{k+1}{k+2} \cdots \frac{N-1}{N}$
$\hspace{3em} = \frac{k}{N}$

データ $S [ i ] \left( k \leq i \right)$ がReservoirに残る確率
$\hspace{3em} = R [ j ]$ に $S [ i ]$ が代入される確率 $\times R [ j ]$ が更新されない確率
$\hspace{3em} = \frac{k}{i} \frac{i}{i+1} \frac{i+1}{i+2} \cdots \frac{N-1}{N}$
$\hspace{3em} = \frac{k}{N}$

実装例

行数が未知のテキストファイルから $k$ 行をReservoir samplingにて取り出す例。
サンプリング結果はテキストファイルとして出力。

import random

def random_sampling(inputFile, k=10000, header=True):
    R = []
    columns = ''
    with open(inputFile) as fin:
        if header: columns = fin.next()
        for i in range(k):
            R.append(fin.next())
        i = k
        for line in fin:
            j = random.randint(0,i)
            if j<k: R[j] = line
            i += 1
    outputFile = inputFile.split('.')[0] + '_random' + str(k) + '.' \
                    + inputFile.split('.')[1]
    with open(outputFile,"w") as fout:
        if header: fout.write(columns)
        for line in R:
            fout.write(line)

参考

Reservoir sampling - Wikipedia, the free encyclopedia

2013-12-25

Twitterと形態素解析

R Advent Calender 25日目の記事を担当しますメロスです。
12月25日である本日は街中がピンク色に包まれるあのクリスマスです。
今年ももうそんな時期なんですね。
ところで、街中がカップルで埋め尽くされる中、
ある重要なプロジェクトが裏ではしっていることをご存知でしょうか。

そう、卒論です。

卒論締め切り間近の学生にとっては、
クリスマスどころではないとか、
いや、むしろ絆を深めるチャンスであるとかないとか。。。

僕は誰かの役に立ちたいという気持ち、
すなわち自己承認欲求がえげつない人間ですから、
彼・彼女らの悲痛な叫びを汲み取って、
こんな僕にも手伝えることはないか考えてみました。

構想としてはこんな感じです。
①Twitterから「卒論」を含むツイートを取得
②ツイートに含まれる形容詞の使用頻度ランキングを出力
③ランキングの結果から卒論に没頭している学生の気持ちを察し、エールを送る

①Twitterから「卒論」を含むツイートを取得

ツイートデータをとってくるためには以下の３つが必要になります。
(1)Twtitterのアカウント
(2)TwitteRパッケージのインストール
(3)ROAuthパッケージのインストール

TwitteRパッケージを利用すると、
ツイートの取得やRからのツイートが可能になります。
詳細はこちらにありますので、ご参考に。
http://cran.r-project.org/web/packages/twitteR/twitteR.pdf

また、ROAuthパッケージはOAuth認証のためのパッケージとなっており、
そのー、深く考えたら負けです。

この記事では(2)以降について記していきます。
なんで、まずはTwitteRパッケージとROAuthパッケージをインストールしましょう。

install.packages("twitteR")
install.packages("ROAuth")

インストールが完了したので、
実際に「卒論」というキーワードを含んでいるツイートを
とってきます。

#「卒論」を含むツイートを取得
library(twitteR)
library(ROAuth)

twit.oauth <- OAuthFactory$new(
	handshakeComplete = TRUE,
	signMethod = "HMAC",
	consumerKey = "*********************",
	consumerSecret = "*********************",
	oauthKey = "*********************",
	oauthSecret = "*********************)
	
registerTwitterOAuth(twit.oauth)

n <- 3200
tweet <- searchTwitter("卒論",n)

これで3200ツイートのデータがとれました。

②ツイートに含まれる形容詞の使用頻度ランキングを出力

ツイートがどういう単語から構成されているのか調べるために
RMeCabを使って形態素解析を行います。
そのためには、以下の２つが必要になります。
(1)MeCabのインストール
　http://mecab.googlecode.com/svn/trunk/mecab/doc/index.html#download
(2)RMeCabパッケージのインストール
　http://rmecab.jp/wiki/index.php?RMeCab

実際にRMeCabを使って、大好きな言葉を形態素解析してみます。

> library(RMeCab)
> RMeCabC("お家に帰るまでが遠足です")
[[1]]
  名詞 
"お家" 

[[2]]
助詞 
"に" 

[[3]]
  動詞 
"帰る" 

[[4]]
  助詞 
"まで" 

[[5]]
助詞 
"が" 

[[6]]
  名詞 
"遠足" 

[[7]]
助動詞 
"です"

RMeCabC()は日本語文字列に対して、
形態素解析の結果をリストとして返す関数です。
「お家に帰るまでが遠足です」という一文が
意味をもつ最小の言語単位に分割されたことになります。

先ほど取得したツイートデータにも形態素解析を行い、
形容詞を対象とした出現回数のランキングを算出してみます。
こんどはRMeCabFreq()を使ってみます。
RMeCabFreq()はテキストファイルを形態素解析し、
その頻度をデータフレームとして返す関数です。

tweet.text <- character(n)
for (i in 1:n){
	tweet.text[i] <- tweet[[i]]$text
	tweet.text[i] <- gsub("RT.?@.*","",tweet.text[i],perl=TRUE)       #RT以降の文章を削除
	tweet.text[i] <- gsub("@\\w*","",tweet.text[i],perl=TRUE)         #アカウント名を削除
	tweet.text[i] <- gsub("http:\\/\\/[\\w\\.\\/]*","",tweet.text[i],perl=TRUE)	#URL削除
}
write(tweet.text,file="tweetdata.txt")

word <- RMeCabFreq("tweetdata.txt")
adjective <- word[word$Info1 == "形容詞",]
adjective <- adjective[order(-adjective$Freq),]

形容詞ランキングトップ30はこんなふうになりました。

> head(adjective,30)
         Term  Info1  Info2 Freq
3955     ない 形容詞   自立   85
3919     いい 形容詞   自立   34
3971   やばい 形容詞   自立   31
4018     眠い 形容詞   自立   27
4004     早い 形容詞   自立   20
4041     いい 形容詞 非自立   20
4007   楽しい 形容詞   自立   19
3952   つらい 形容詞   自立   14
3992   忙しい 形容詞   自立   13
4013     無い 形容詞   自立   11
3946   すごい 形容詞   自立   10
3999     悪い 形容詞   自立   10
4027     良い 形容詞   自立   10
3972     よい 形容詞   自立    9
4032     辛い 形容詞   自立    8
3985     多い 形容詞   自立    7
3986   嬉しい 形容詞   自立    7
4015     痛い 形容詞   自立    7
4034     遅い 形容詞   自立    7
3935 かわいい 形容詞   自立    6
3958   はやい 形容詞   自立    6
3984   可愛い 形容詞   自立    6
3993     怖い 形容詞   自立    6
3916   っぽい 形容詞   接尾    5
3937   きつい 形容詞   自立    5
4008   欲しい 形容詞   自立    5
4047     良い 形容詞 非自立    5
3922   うまい 形容詞   自立    4
3928 おかしい 形容詞   自立    4
3949   だるい 形容詞   自立    4

あー、なるほどね。うん。
わかる。わかる。
大丈夫です。
完璧に理解しました。

卒論で苦しんでいる彼・彼女らへのエールでしたら、
僕に任せてください。

③ランキングの結果から卒論に没頭している学生の気持ちを察し、エールを送る

進捗どうですか♡

2013-08-04

ラグランジュの未定常数法（the method of Lagrange mulitipliers)

mathematics

概要

多変数実数値関数の停留点を変数間の制約条件下で求めるための方法。

束縛問題を極値問題として扱うことができる。

そのため、束縛されたいと強く望むあなたには、おすすめできない。

等式制約条件下での最適化問題

目的関数 $\hspace{3} f( \bf{ x }) \hspace{3}$ が最大となる点を

$m$ 個の等式制約

$h_i(\bf{x})=0,\hspace{10} i=1,2,\cdots,m$

の下で求める。

この問題はラグランジュの未定常数法により解くことができる。

※等式制約条件下での最適化問題をラグランジュの未定常数法で解く場合、目的関数や制約条件式に対して微分可能性や凸関数である等の条件が必要であるが、それらがすべて満たされていると仮定する。

この問題に対してラグランジュ関数と呼ばれる次の実数値関数を定義する。

$L(\bf{x},\bf{\lambda}) = f(\bf{x}) + \displaystyle\sum_{i=1}^m \lambda_ih_i(\bf{x})$

このとき、以下の連立方程式を解くことによって最適点 $\hspace{3}\bf{x}_*\hspace{3}$ を求めることができる。

$\frac{\partial L(\bf{x},\bf{\lambda})}{\partial \bf{x}} = \frac{\partial f(\bf{x})}{\partial \bf{x}} + \displaystyle\sum_{i=1}^m \lambda_i\frac{\partial h_i(\bf{x})}{\partial \bf{x}} = \bf{\0}$

$\frac{\partial L(\bf{x},\bf{\lambda})}{\partial \lambda_i} = h_i(\bf{x}) = 0, \hspace{5} i=1,2,\cdots,m$

このようにして最適点を求める魔法が昨今、世間を賑わせているラグランジュの未定常数法である。

幾何学的な説明

例題を通して、幾何学的な説明をしていく。

例：

$max.\hspace{10} f(x_1,x_2) = -x_1^2-x_2^2+8$

$s.t.\hspace{10} h(x_1,x_2)=x_1+x_2-8=0$

このとき、目的関数 $\hspace{3}f\hspace{3}$ のグラフは以下のようになる。

何かに似ているが偶然である。

f:id:dr_melos:20130804180546p:plain

これを上から見た図が以下。

※黒い線は等高線を、青い線は制約によって $\hspace{3}(x_1,x_2)\hspace{3}$ がとりうる範囲（制約面）をそれぞれ表している。

f:id:dr_melos:20130804190552p:plain

まず、最適点 $\hspace{3}\bf{x}_*\hspace{3}$ においては制約面と等高線は接していなければならない。

もしそうでないとすると、制約面に沿って $\hspace{3}f(\bf{x})\hspace{3}$ の値がさらに大きくなるように移動できてしまい、矛盾。

よって以下のことがいえる。

$\exists\lambda \neq 0,\hspace{10}\nabla_{\bf{x}}f(\bf{x}) + \lambda\nabla_{\bf{x}}h(\bf{x})=0$

また、制約面と等高線が一致していたとしても、制約面上に解があるとは限らないので、

$\begin{equation} h(\bf{x})=0\end{equation}$

という条件も必要になる。

以上をまとめると次の連立方程式が導かれる。

$\begin{eqnarray} \left\{ \begin{array}{l} -2x_1 + \lambda=0\\ -2x_2 + \lambda=0\\ x_1+x_2-8=0 \end{array} \right. \end{eqnarray}$

これより、最適解 $\hspace{3} $x_1,x_2$ = $2,2$ \hspace{3}$ が得られた。

不等式制約条件下での最適化問題

ラグランジュの未定常数法は不等式制約条件下での最適化問題にも適用できる。

ただ、僕はもう疲れ（手記はここで途絶えている）

参考

パターン認識と機械学習上

言語処理のための機械学習入門

2013-07-14

MySQLをMacOSにインストール

MySQL

インストール

①ここからファイルをダウンロード

②ダウンロードしたファイルを開いて
　次の２つのインストーラパッケージを実行
・mysql-5.6.12-osx10.7-x86_64.pkg
・MySQLStartupItem.pkg

③MySQL.prefpane　を実行して
　システム環境設定にMySQLを追加

これでインストールは終わり

MySQLの起動と終了

起動時にはターミナルから次のコマンドを入力

sudo /Library/StartupItems/MySQLCOM/MySQLCOM start
mysql -u root -p

パスワードを聞かれたらEnter

終了時にはターミナルから次のコマンドを入力

sudo /Library/StartupItems/MySQLCOM/MySQLCOM stop

2013-01-29

Rで行列を使って累積を求める方法

Rでベクトルの累積和を出す際に
行列かけてできないかな
なんてことを考えてみました。

#もしもこんなベクトルがあったとしたら
x <- c(1:10)

#ベクトルの大きさを調べて、その大きさの正方行列を作成
n <- length(x)
A <- matrix(1,n,n)

#上三角行列を0にする
A[upper.tri(A)] <- 0

#Aはこんな行列
> A
      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
 [1,]    1    0    0    0    0    0    0    0    0     0
 [2,]    1    1    0    0    0    0    0    0    0     0
 [3,]    1    1    1    0    0    0    0    0    0     0
 [4,]    1    1    1    1    0    0    0    0    0     0
 [5,]    1    1    1    1    1    0    0    0    0     0
 [6,]    1    1    1    1    1    1    0    0    0     0
 [7,]    1    1    1    1    1    1    1    0    0     0
 [8,]    1    1    1    1    1    1    1    1    0     0
 [9,]    1    1    1    1    1    1    1    1    1     0
[10,]    1    1    1    1    1    1    1    1    1     1

#各成分までの累積を求める
y <- A %*% x

#結果はこんな
> y
      [,1]
 [1,]    1
 [2,]    3
 [3,]    6
 [4,]   10
 [5,]   15
 [6,]   21
 [7,]   28
 [8,]   36
 [9,]   45
[10,]   55

2012-12-09

序章：どうも僕です。

そうです。僕があのメロスです。

昔は約束を守るために走ったりしましたが、

もうね、

この年になると自分が一番大事ですよね。

走ってらんないわけです。

セリヌンティウスとかどうでもいいし。

これが大人になるということなんでしょうか。

好きなことをやって生きていけたら

それは素晴らしいことで、

それを実現するために

流した汗（あせ）とか汚い汁（しる）とかを

一度、浄水器を通してから

ここに残していこうかなと企んでいます。