TIM Labs
このところ、とても忙しくなって、書くペースが落ちてしまった。
というのも、R&D部門が新しいオフィスに移動作業中だからなのだ。
真新しいビルが用意され、そのビルの一室に入ることになった。

場所は、なんと、電気通信大学の100周年記念キャンパスの中になった。
キャンパス名に100周年記念キャンパスということで、どこか東京の郊外、不便なところにキャンパスが出来たわけではなく、元のキャンパスと道一つ(甲州街道)隔てた場所にある。

電気通信大学(電通大)については、いまさら説明することもないだろう。 Wikipediaで調べるのも良いと思う。
アンリツの前進である安中電機製作所内に1916年にできた無線講習所が本来の起源だと思うが、公的機関としてスタートしてから100年を100周年ということになっている。
無線が一気に世界の注目を集めたのはタイタニック号の沈没により、無線の重要性、情報通信の重要性が認識されたのだった。
それ以来、無線&通信はどんどん進歩し、今のインターネット時代を支えている技術である。
有線の通信だけだったら、インターネットはここまで幅広く使われることはなかっただろう。
歴史的なことは、電気通信大学の同窓会である目黒会ウェブマガジンの記事が詳しい。

ということで、ITの根幹の技術をずっと担ってきた電通大の中に研究拠点を設けて、新しいことができるかなと思っているところだ。

DSCN0481 (450x600).jpg DSCN0473 (600x450).jpg
キャンパス、建物はこんな感じで、ずいぶん立派なんだ。
2017年春に竣工し、この春から使い出したので、ピカピカである。

新しいだけではなく、実は、色々面白い仕掛けというか、いろいろあるのだが、今回は1つだけ紹介しよう。

ネット上でも話題になっている通り、ここの1階には電子パーツを売る24時間営業のコンビニ・セブンイレブンが入っている。

DSCN0480 (600x450).jpg DSCN0267 (450x600).jpg
さすがに電通大内のコンビニだけあって、工具やテスター、パーツなどが置かれているコーナーがあるのだ。
調布、府中、、、などで夜に電子部品が欲しいとなったら、秋葉原まで行かなくても、調布で入手可能なのだ。
甲州街道に面しているので、とても分かりやすいと思う。

実は、これだけではなかったのだが、それについては別の機会に紹介しよう。


Pythonには、様々なモジュールがある。
それを1つずつ毎週紹介すると何年も話が持ちそうだが、こちらの気力が持ちそうにない。実際、それほどたくさんのモジュールが存在するのだ。

Excelをいじるためのモジュールもいくつかあるようなのだが、今回はopenpyxlを紹介する。

openpyxl Documentationという立派なドキュメントがある。
英語のドキュメントで、311ページもあるので、真面目に読むのは大変だ。

A Python library to read/write Excel 2010 xlsx/xlsm filesがもう少し楽そうなドキュメントである。

実際には、東邦大学理学部のPythonからExcelファイルをいじるopenpyxl
あたりを出発点としてあれこれ探して次のプログラムをでっち上げただけである。


#!/usr/bin/env python

from PIL import Image, ImageColor
import openpyxl
from openpyxl.styles import PatternFill
from openpyxl.utils import get_column_letter

def coltostr(coltuple):
    r, g, b = coltuple
    a = 255
    return format(((a*256+r)*256+g)*256+b, '08X')

def image2excel(imagefile,excelfile,sheetname):
    im = Image.open(imagefile)
    width, height = im.size
    
    wb = openpyxl.load_workbook(excelfile)
    sheet = wb.get_sheet_by_name(sheetname)
    
    # セルサイズの調整
    for x in range(0,width):
        sheet.column_dimensions[get_column_letter(x+1)].width = 1
    for y in range(0,height):
        sheet.row_dimensions[y+1].height = 5.7
    
    # 画像ファイル ⇒ Excelワークシート
    for x in range(0,width):
        for y in range(0,height):
            colstr = coltostr(im.getpixel((x,y)))
            cell = sheet.cell(row=y+1, column=x+1)
            cell.fill = PatternFill(fill_type='solid',fgColor=colstr)
    
    wb.save(excelfile)
    im.close()

image2excel("アンパン.jpg","アンパン.xlsx","アンパン")
暑くなってきて頭が回らない。
ということで、カフェに逃げ込んで、こんなものを食べながら、どうでも良いことを考えた。

ampan300.pngすると、ついこんなイメージが頭に浮かんでしまった。

アンパンonExcel.png
画像の周りにいっぱい罫線ができている。
つまり、これは、Excel、OpenOfficeのCalcとか、その類だ。

ということで、Excelの各セルに、元のイメージファイルの画素を対応させてみよう。

アンパンの画像は、元は300x300だったのだが、いかにもセルに色を付けた感が出るように、64x64とそれなりに粗く(小さく)した。

アンパン.jpgさて、どうやって描こうか。
Excelだから当然VBAのプログラムを組めばよい、と思いたいが、VBAのプログラムは作りたくない。
どうしよう。

このところ、ずっとPythonをいじっているのだが、星の数ほどあるPythonのモジュールの中には、Excelのモジュールも当然あるはずだ、きっとある。

ということで、次回までに作ってみよう。
なお、アンパンは美味しかった。


NLTK (313x400).jpg

題名:入門自然言語処理

Natural Language Processing with Python
Analyzing Text with the Natural Language Toolkit

Steven Bird、Ewan Klein、Edward Loper 著、萩原 正人、中山 敬広、水野 貴明 訳
大型本、592項、本体3800円

2010年11月8日 発行

オライリー・ジャパン

Natural Language Toolkit

Natural Language Processing with Python (書籍)


人工知能、機械学習、ディープラーニングというと、画像処理関連がやたらに多いが、それ以外の分野もある。
その中でも、自然言語処理は非常に大きな、そして重要な分野である。

ことばをコンピュータで扱おうとすると、画像とは違った、あれこれ面倒なことがいっぱいある。
言葉を対象としている人工知能の本の場合、自然言語処理の部分の説明は非常に短く、いきなり読もうとしても用語が分からない、どんなツールがあるのか、サンプルデータがあるのか、だいたい分からないことだらけになる。

自然言語処理を対象としてAI関連の本で、自著で延々と説明するのは大変なので、読むべき自然言語処理の本が挙げられていることが多いが、そのなかで必ずといってよいくらい紹介されるのが、この『入門自然言語処理』である。

この本は、Pythonを使って、自然言語処理の基本を紹介している。
といっても、原書は英語で、日本語の場合どうなんだろうと思ったら、最終章が「Pythonによる日本語自然言語処理」となっている。

さて、この本、発行が2010年とかなり古く、原書は、2009年になっている。
そのため、Python3ではない。

この本で使われているのが、NLTK(Natural Language Toolkit)という、Pythonのツールキットである。
このツールキットは、アメリカのアイビーリーグの1つ、ペンシルベニア大学にて作られたものだ。


以上は前置きで、これから肝心なことを紹介しよう。

本書はとても古いのだが、Natural Language Toolkit のサイトでは、今も更新が続いており、ちゃんとPython 3 対応になっている。

ソフトだけでなく、書籍の方も、ネット上はちょこちょこと更新されているように見える。
さらに、これらは、オープンであり、自由に使えるので、とても助かる。

本は、文字だけでなく、プログラムや実行例が多数載っており、これらを自身のPython上で確かめるには、オンラインの書籍からコピペをいっぱいすることで、確認ペースも上昇する。

オンライン版は文章の部分は英語であるが、肝心なのはプログラム、実行例などであろう。
それほど大きくは変わっていないようなので、英語をどうしても読みたくない場合には、英語版を見ながら、文章を読むときだけ翻訳書に頼るという方法もある。

でも、結局面倒になるので、全部オンラインだけで済ませるのが効率がよく、かつお金もかからない。

英語を勉強するのではなく、英語で勉強しよう。
金もかからず、技術も身に付き、情報はいっぱい集まる。


前回、SATについて簡単過ぎる紹介をした。
もっと詳しく知りたい場合は、神戸大学情報基盤センターの田村教授の資料 情報基礎特論:制約プログラミング入門 (PDF, 98ページ)が相当詳しい。

SATには、制御文がない。つまり、if-else, while, for,...などの文が存在しない。
存在するのは、定義文、条件文だけである。
つまり、定義の範囲で、条件を全て満たす変数の値を探すだけである。
もちろん、条件を満たす解が1組以上ある場合もあるが、解の個数を、解なし(0)、解1(ユニーク解)、複数解(多重解)の区別もしてくれたりする。

さて、こんなSATは、どうやって問題を解いているだろうか。
結論から言うと、最終的には虱潰しをやっているだけである。
しかし、途中で検索を省略できる箇所、早々に決まってしまう箇所、選択肢が非常に狭められる場合などを考慮し、相当高速に虱潰しができるようになっている。

それでも、さすがに25x25のナンプレを調べつくすのには相当の時間がかかると予想した。

神戸大学の田村教授に16x16の問題を解いてもらったときの時間が
real    0m3.311s
user    0m5.244s
sys     0m0.280s
20x20の問題を解いてもらったときの時間が
real    0m5.079s
user    0m7.312s
sys     0m0.600s
であった。
組み合わせ問題なので、サイズが大きくなると急激に時間がかかるようになるものだが、マスの総数の増加割合よりやや時間がかかる程度であった。 これなら、25x25も解かれてしまうかも知れない。
25x25のヒント数の少ない問題を作るのは、非常に困難で、コンピュータリソースをやたらに食う。

25x25の問題を送ったら、1時間余したら返信があり、なんと
real    0m17.693s
user    0m21.312s
sys     0m0.332s
という短い時間の間に解かれてしまったのだった。 最近のSATの性能向上には目を見張るものがある。 SATは、明示的な、手筋のような解き方をプログラミングする必要がないので、解法アルゴリズムがまったく思いつかないような場合に有効な方法だ。 このまま高性能化が進むと、アルゴリズムを考えず、とりあえずSATに探させれば良い、と考えるようになりそうだ。

いや、SATの開発者は、効率向上のために、ものすごくアルゴリズムを最適化しているはずだ。

ナンプレ(あるいはSUDOKU)というパズルをご存知だろうか。
通常は、9x9の問題なのだが、つい、25x25の問題を作った者がいる。
問題に最初から示されている数のことをヒントと言い、その個数をヒント数と言うのだが、以下がヒント数146(現在世界最少)の問題である。

+--------------+--------------+--------------+--------------+--------------+
| . 12  .  .  .| .  .  .  .  .| .  .  .  9  .| .  . 15  .  .|22  .  .  .  .|
| .  .  .  .  .| .  9  . 19  .| .  . 10 11  .| .  .  .  .  .| .  .  .  .  .|
| .  4  . 22  .| .  .  .  .  .| .  .  .  .  .| .  . 12  .  .|20 15  1  .  .|
|16  1 20 15  .| .  .  .  .  .| .  .  .  .  .|14  .  4  . 22|12 25  .  .  .|
| .  .  .  .  .| .  7  2 11  .|23  . 19  8  .| .  .  . 13  .| .  .  .  .  .|
+--------------+--------------+--------------+--------------+--------------+
|13  .  8  .  2| .  .  .  .  .| .  .  7 23  6| .  9  . 19 11| .  .  .  .  .|
| .  .  .  . 23| .  .  .  . 16| .  .  .  .  .| .  .  .  .  .| 1  .  .  .  .|
| 7  .  .  . 10| 3  .  .  .  .| .  .  9 19  .| . 13  . 23  .| .  .  .  5  .|
| .  .  .  .  .|15  .  .  . 22| .  .  .  .  .| .  .  .  .  .|25 20  .  .  .|
| .  .  .  .  .|12  . 14  1 25| .  .  .  .  .| .  .  3  .  .|16  4 15  .  .|
+--------------+--------------+--------------+--------------+--------------+
| .  .  .  .  .| . 19  9  .  .| .  . 13  7  .| .  .  .  5  .| .  .  . 23 10|
| . 22  . 25 17| .  .  .  .  .| .  .  .  .  .|12  . 20  .  .| .  .  .  .  .|
| . 20 12 16  .| .  .  .  .  .| .  .  .  . 14|15 22  1  . 25| .  .  .  .  .|
| . 15  .  .  .| . 11  .  .  .| .  .  .  .  .| .  . 16  .  .| .  .  .  9  .|
| .  .  .  1  .| . 10  . 23  .| .  .  .  . 18| .  .  .  .  .| .  .  .  .  8|
+--------------+--------------+--------------+--------------+--------------+
|10  .  .  .  8| . 13  .  5  .| .  .  .  .  .| . 19  . 11 23| .  .  .  6  .|
| .  .  . 17  7| .  .  .  .  .| .  .  .  .  1| .  .  .  .  .| 4 22  .  .  .|
| .  .  .  . 11| . 23  .  .  .| .  .  .  . 20| .  .  .  2  .|14  .  .  .  .|
|19  . 23  .  5| .  8  .  9  .| . 21  .  .  .| . 10  .  7  .| .  .  .  .  .|
| .  3  .  .  .| .  .  .  .  .|25  4  .  . 12| .  .  .  .  .|15  1 16  .  .|
+--------------+--------------+--------------+--------------+--------------+
| .  .  .  .  .| .  .  .  . 15| . 12  .  . 25| 1  . 22  .  .| 3  .  .  .  .|
|23  .  .  . 19| .  2  .  .  .| .  .  .  .  .| .  .  . 10  .| .  .  .  7 11|
| .  .  . 18  .| .  .  .  .  .| . 20  .  .  .| .  .  .  .  .| .  .  .  .  .|
| .  .  .  .  .| .  .  .  .  4|14 15  .  . 22| .  .  .  .  .| .  .  . 10  .|
|11  .  .  .  9| .  .  .  .  .| .  .  .  .  .| .  .  .  .  .| .  .  . 19  .|
+--------------+--------------+--------------+--------------+--------------+
ナンプレの問題を解くプログラムはネット上にゴロゴロ転がっているのだが、大きいサイズに対応したソルバーはなかなか存在しない。
 小さいときには、手筋を実装せず、単に再帰だけでも十分高速に解けるのだが、それではサイズが大きくなると直ぐに破綻するはずだ。

この問題、なかなか解いてくれる人がいない。
ナンプレに限らないのだが、パズルの多くは、限られたルールとヒントから解くのである。
そして、こういう問題のことを、制約充足問題と言う。
日本語での説明は貧弱なので、ぜひ Constraint satisfaction problem を参照されたい。

さて、制約充足問題というと、すぐに思いつくのがSATであろう。
ということで調べると、SATでパズルを解く研究をしている神戸大学情報基盤センターが直ぐに見つかる。

ということで、上の25x25の問題を送ってみた。さて、どうなるだろうか?

このブログ記事は『データ解析のための統計モデリング入門』(久保拓弥 著、岩波書店)という、とても分かりやすい統計モデリングの入門書を、さらに分かりやすくするために読書メモをまとめたものです。 今回は第8章、マルコフ連鎖モンテカルロ(MCMC)法とベイズ統計モデルについてのまとめの四回目です。

この章ではMCMCやメトロポリス法について説明がされています。 その中で、メトロポリス法によるサンプリングと尤度が関係あることが解説されています。 なので、メトロポリス法によるサンプリングの分布と尤度を比べられるコードを用意しました。 コードはRで書きました。

dataSize <- 20
N <- 8

data <- rbinom(dataSize, N, 0.45)

likelihood <- function (q) {
  sum(log(dbinom(data, N, q)))
}

develop <- function(q) {
  if (runif(1, 0, 1) < 0.5) {
    q + 0.01
  } else {
    q - 0.01
  }
}

sampling <- function(n, init) {
  q <- init
  ql <- likelihood(q)
  samples <- vector(length = n)
  i <- 1
  while(i <= n) {
    p <- develop(q)   
    pl <- likelihood(p)
    r <- log(runif(1, 0, 1))
    if (ql < pl || r < pl - ql) {
      q <- p
      ql <- pl
    }
    samples[i] <- p
    i <- i + 1
  }
  samples
}

qs <- seq(0, 1, by = 0.01)
ls <- vector(length = length(qs))
for(i in 1:length(qs)) {
  ls[i] <- prod(dbinom(data, N, qs[i]))
}
plot(qs, ls, type = "l", xlab = "q", ylab = "Likelihood")

samples <- sampling(10000, 0.3)
hist(samples, breaks = seq(0, 1, by = 0.01), main = "10000 samples", xlab = "q")

サンプルデータとしては、本の中でやっているのと同じように、事象の発生確率が 0.45 で、事象の試行回数は 8 の二項分布から、20 個のサンプルを取って用意しました。 コードを実行すると次のような図がプロットされます。

metro_likelihood.png

この図は q を変化させた時の尤度のグラフです。 発生させたサンプルデータの偏り方によって多少ずれますが、だいたい 0.45 を中心とした山なりなグラフになります。

metro_hist.png

この図はサンプル数 10000 の場合のメトロポリス法によるサンプリングの結果のヒストグラムです。 確かに、さきほどプロットした尤度に比例した分布が得られていることが分かります。

ダメな統計学 (283x400).jpg

題名:ダメな統計学

悲惨なほど完全なる手引書

アレックス・ラインハート 著  西原 史暁 訳

A5、185項、本体2200円

2017年1月20日 発行

勁草書房

https://www.statisticsdonewrong.com/ (オリジナルサイト)

http://id.fnshr.info/2017/01/20/sdw-trans-publ/ (日本語サイト)


最近、世の中、統計学が流行っているようである。

確かに昔と違って、大量のデータもさくっと処理できるだけの計算機パワーがあり、様々な分析、予想、判断などに実際に利用できるようになった。

そして、ビッグデータ、人工知能などのためには、統計学は当然知っていなければいけない。


ということで、統計学を勉強と思ってではなく、世の中、どう考えてもオカシイ、誤った統計の利用、誤用、乱用、悪用が目に付いていたところでこの本が出版されたので読んでみた。


第一印象は、結構読みづらい本であった。

日本語訳があまりよくない。元々内容がねちねちしたところがあるのは当然なのだが、それゆえしっかりした、かつ読みやすい日本語にして欲しいところだが、その当たりはかなり改善して欲しいところだ。英語がすらすら読めるのならば、英語で読むべきかもしれない。

内容は、論文などでの統計の利用がいかに正しくないかが書かれている。
世界の主要な論文掲載誌に載る論文も、1/3程度は統計的におかしなことを書いており、結論が違っているとか。

とくに、医学、医療関連、つまり治療の有効性の判定などでの誤りの指摘が多い。
医療系の世界では、統計の誤解だけでなく、どうやら政治、ビジネス的な問題がからんで、世界中どこもドロドロなところがあると。

統計とくに検定についての言及が多いが、そのあたりを正しく把握しようと思うと、一通りの統計学の知識が必要になる。

いろいろ統計的に誤った実例が載っている。とくに、薬で、有効と認定されていたのが無効と判定されてしまったいきさつなどの紹介もある。

こういう研究は、直接的な成果に結びつかないので、関心を持たれることがないが、統計にきちんと基づいて論文の査読がされれば、誤った結論を出してしまうことが相当減少することは確実だ。

一部の論文誌では、査読に統計、検定の根拠を示すのを義務化するのも進み始めたみたいだ。

本書は、さくっと読んで、しっかり統計学を勉強するのが重要だ。


なお、本書の内容であるが、英語も、翻訳もネットでフリーで読むことができる。

このブログ記事は『データ解析のための統計モデリング入門』(久保拓弥 著、岩波書店)という、とても分かりやすい統計モデリングの入門書を、さらに分かりやすくするために読書メモをまとめたものです。 今回は第8章、マルコフ連鎖モンテカルロ(MCMC)法とベイズ統計モデルについてのまとめの三回目です。

この章ではMCMCやメトロポリス法について説明がされています。 その中で、メトロポリス法によるサンプリングが解説されています。 なので、同じようにメトロポリス法でサンプリングしてみるコードを用意しました。 コードはRで書きました。

dataSize <- 20
N <- 8

data <- rbinom(dataSize, N, 0.45)

likelihood <- function (q) {
  sum(log(dbinom(data, N, q)))
}

develop <- function(q) {
  if (runif(1, 0, 1) < 0.5) {
    q + 0.01
  } else {
    q - 0.01
  }
}

sampling <- function(n, init) {
  q <- init
  ql <- likelihood(q)
  samples <- vector(length = n)
  i <- 1
  while(i <= n) {
    p <- develop(q)   
    pl <- likelihood(p)
    r <- log(runif(1, 0, 1))
    if (ql < pl || r < pl - ql) {
      q <- p
      ql <- pl
    }
    samples[i] <- p
    i <- i + 1
  }
  samples
}

samples1 <- sampling(500, 0.3)
hist(samples1, breaks = seq(0, 1, by = 0.01), main = "500 samples", xlab = "q")

samples2 <- c(samples1, sampling(500, samples1[500]))
hist(samples2, breaks = seq(0, 1, by = 0.01), main = "1000 samples", xlab = "q")

samples3 <- c(samples2, sampling(1000, samples2[1000]))
hist(samples3, breaks = seq(0, 1, by = 0.01), main = "2000 samples", xlab = "q")

samples4 <- c(samples3, sampling(8000, samples3[2000]))
hist(samples4, breaks = seq(0, 1, by = 0.01), main = "10000 samples", xlab = "q")

サンプルデータとしては、本の中でやっているのと同じように、事象の発生確率が 0.45 で、事象の試行回数は 8 の二項分布から、20 個のサンプルを取って用意しました。 なので、このデータに対して q のサンプリングをすると、だいたい 0.45 を中心とした分布が得られるはずです。

コードではサンプル数が 5001000200010000 の時にヒストグラムをプロットしています。 コードを実行すると次のような図がプロットされます。

metro_500.png

この図はサンプル数 500 の場合です。 まだヒストグラムが凸凹しているのが分かります。

metro_10000.png

この図はサンプル数 10000 の場合です。 だいぶヒストグラムがなめらかになって、分布が収束してきているのが分かります。

PythonをAIなどで使おうと思ったら、かならず使うことになるのがNumPyである。 様々な計算処理をしようとすると、はやり配列を使って高速に計算する必要がでてくるので、配列をサポートしているNumPyは避けて通れない。
 NumPyにあるtutorilaなどを参考に勉強したり、ネット上にウヨウヨ存在する説明をみて勉強するのもあるのだが、もうちょっと面白い勉強方法、あるいは腕試しには、100 numpy exercises が良いかも。

題名どおり、100の小問で構成されている。
難易度が3段階になっていて、やさしい、基本的なものから、次第に難しい問題になっているようだ。
最初は1、2行で十分な簡単な問題ばかりだが、次第に長いプログラムを書かないといけなくなり、20行程度のプログラムになることもある。

問題と解答があり、どうしてそう書けば良いのかは自分で考えよう、というスタンスのようである。
もちろん、世の中には延々と解いて、解説ページを公開している人もいる。

解説ページを読むより、解答が理解できなかったら、ネットで色々探して、理解できるようになろう。
さらに、出題範囲以外のことを色々試すことをお薦めする。

あ、書き忘れていたことがある。
もちろん、全部英語だけれど、やさしいプログラミング英語なので、問題ないだろう。
どうせ、プログラムをいっぱいやるようになったら、英語しかない世界にすぐなるのだし、プログラミング英語に慣れる良い機会だ。

まず、最初の問題。こんなに簡単で、分からない人はいないはずだ。

1. Import the numpy package under the name np (★☆☆)

import numpy as np

そして最終問題は、こんな感じだ。

100. Compute bootstrapped 95% confidence intervals for the mean of a 1D array X (i.e., resample the elements of an array with replacement N times, compute the mean of each sample, and then compute percentiles over the means). (★★★)

# Author: Jessica B. Hamrick

X = np.random.randn(100) # random 1D array N = 1000 # number of bootstrap samples idx = np.random.randint(0, X.size, (N, X.size)) means = X[idx].mean(axis=1) confint = np.percentile(means, [2.5, 97.5]) print(confint)


こんな問題に対するプログラムをすらすら書けるようになったら、もちろんNumPyのお勉強は優秀な成績で卒業だ。
NumPyの勉強は、通過すべき最初のゲートみたいなものだ。どんどん先へ行こう。

最近のコメント