ほくそ笑む

R言語と統計解析について

「子供に解けて大人に解けない問題」を統計的に無理やり解いてみた

今日は、R-bloggers に面白い記事が上がっていたので、それを紹介してみようと思います。

問題

「子供にはすぐに解けて、大人にはなかなか解けない不思議な問題」をご存知でしょうか?
最近ネットで割と話題になりました。
その問題は、次のようなものです。

8809 = 6
7111 = 0
2172 = 0
6666 = 4
1111 = 0
3213 = 0
7662 = 2
9312 = 1
0000 = 4
2222 = 0
3333 = 0
5555 = 0
8193 = 3
8096 = 5
7777 = 0
9999 = 4
7756 = 1
6855 = 3
9881 = 5
5531 = 0
2581 = ?

https://twitter.com/#!/yappyJP/statuses/172086299099004928

なかなか面白い問題です。
答えはここでは書きませんので十分楽しんでから下を読んでください。
本題に入りますが、この問題を統計的に解こうと試みたブログ記事がありました。
つまり、問題で例示されたデータをもとに、未知の入力に対する解答を統計的に算出しようというのです。
Solving easy problems the hard way | R-bloggers
この記事は結果もなかなかおもしろかったので、今日はこれを紹介しようと思います。

データを変換する

この問題、このままだと統計手法にうまく当てはめることができないので、左辺にある 0-9 の数をカウントしたデータに変換します。
例えば、次のような感じです。

8809 = 6 -(変換)-> 1 0 0 0 0 0 0 0 2 1 = 6
7111 = 0 -(変換)-> 0 3 0 0 0 0 0 1 0 0 = 0
2172 = 0 -(変換)-> 0 1 2 0 0 0 0 1 0 0 = 0

ちょっと分かりにくいかもしれないので、実際にデータ変換プログラムを統計言語 R で書いてみましょう。
もともとのデータはここにあります。

series <- read.csv("https://raw.github.com/gist/2061284/44a4dc9b304249e7ab3add86bc245b6be64d2cdd/problem.csv")
freqTable <- as.data.frame( t(apply(series[,1:4], 1, function(X) table(c(X, 0:9))-1)) )
names(freqTable) <- c("zero","one","two","three","four","five","six","seven","eight","nine")
freqTable$dep <- series[,5]
head(freqTable)

変換後のデータの最初の6行だけ見てみましょう。

  zero one two three four five six seven eight nine dep
1    1   0   0     0    0    0   0     0     2    1   6
2    0   3   0     0    0    0   0     1     0    0   0
3    0   1   2     0    0    0   0     1     0    0   0
4    0   0   0     0    0    0   4     0     0    0   4
5    0   4   0     0    0    0   0     0     0    0   0
6    0   1   1     2    0    0   0     0     0    0   0

8809 は、0 が 1個、8 が 2個、9 が 1個、その他は 0個なので、1 0 0 0 0 0 0 0 2 1 となります。
dep は右辺そのままです。
7111 = 0, 2172 = 0, 6666 = 4, 1111 = 0, 3213 = 0 も同様に変換されています。
このように変換したデータを統計解析にかけてみようと思います。

線形回帰モデル

この問題を、具体的には線形回帰モデルで解いてみます。
線形回帰モデルというのは、答えが入力変数の線形結合で表わされるという仮定です。
例えば、3つの変数 x1, x2, x3 があった場合、答え y は

y = a1 * x1 + a2 * x2 + a3 * x3

として表されます。
ここで、a1 〜 a3 は係数で、これを求めることが問題を解くことにあたります。
今、変数は zero から nine まで、10個あるので、答え dep は次のように表されます。

dep = a0 * zero + a1 * one + a2 * two + ...... + a9 * nine

このときの係数 a0 〜 a9 を最小二乗法を用いて計算してみます。

myModel <- lm(dep ~ 0 + zero + one + two + three + four + five + six + seven + eight + nine, data=freqTable)
round(myModel$coefficients)

結果は以下のようになりました。

 zero   one   two three  four  five   six seven eight  nine 
    1     0     0     0    NA     0     1     0     2     1 

four だけ問題文に出てこないので、NA となりましたが、あとの係数は求まっています。
これをさっきのモデルの係数に当てはめると、

dep = 1*zero + 0*two + 0*three + 0*five + 1*six + 0*seven + 2*eight + 1*nine

0 の項を省略すると、

dep = zero + six + 2*eight + nine

となります。
日本語で書けば、

答え = 0の数 + 6の数 + 8の数×2 + 9の数

ですね。
この結果に現れているのは内部に○を持つ数だけであり、8については2倍されています。
答えをご存じの方にはピンとくるはずです。
つまり、この結果は答えが数字内部の○の数を加算したものであるということを示しています。
完璧な解答が得られました!

おわりに

ちなみに、問題文の最後の 2581 は 0 + 0 + 1*2 + 0 = 2 という風に求めることができます。
もちろん他の数字に関しても、この式に当てはめれば答えを求めることができます。
こういうクイズを統計手法を使って解くというのは面白い試みですね。
それではまた。