ラベル 散布図 の投稿を表示しています。 すべての投稿を表示
ラベル 散布図 の投稿を表示しています。 すべての投稿を表示

2012/09/18

文系のための「二変数の関係」(2)

データを観察する上で最も重要なことは、何度でも言うけれども、
データの中心」と「データのバラツキ」を把握すること。
あらゆる統計的な手法はこの二つの点が全ての原点になっている。

そう言えば、数量化理論を創り上げた林知己夫氏が、
以下のように語っていることを思い出した。

科学は所詮平均値の議論である。・・・(中略)・・・
平均や構造をみることで大きな知見を得ることができ、
さらに個々を介して新たに発展の道をたどることになる。

この言は、日本行動計量学会35年記念誌『行動計量学』2008年9月号の、
「論文 ーその世俗的な話しー」と題したエッセイに書かれてある。
思わず、苦笑いしたくなる名言が詰まっていて、是非、研究者には読んでもらいたい。

初っ端から話が逸れてしまったが、重要なことは、あくまでデータの中心であって、
その中心からのバラツキを測った指標として、分散標準偏差があった。

さらに、このバラツキ二つの変数間で観察するための指標として、
共分散というものがあり、それを基準化したものとして相関係数があった。
ふむふむ。では、相関係数とはどのようなものだったか?

確か、二つの変数の関係の強さ-1〜1の範囲で表し、

一方が増加するともう一方が減少する場合を負の相関と呼び、
一方が増加するともう一方も増加する場合を正の相関と呼んだ。

そうそう。そういった関連性の強さを表すものであった。

実を言うと、このように二つの変数の関連の強さを観察する方法というのは、
相関係数だけではない。別の考え方に基づく方法がある。
今回の話は、その方法についての話。

とりあえず、相関の話で使ったデータの読み込みから。

# Windows と Linux の人は次のコマンド
X <- read.table("clipboard", sep=",", header=TRUE)

# Mac の人は次のコマンド
X <- read.table(pipe("pbpaste") , sep=",", header=TRUE)

以下が、そのデータ。前回のものを再掲しただけ。
データの説明に関しては、平均の話を参照。

Oroshi_Kg, Taka, Naka, Yasu
Tai, 1336, 3150, 1217, 53
Chinu, 226, 630, 513, 105
Sawara, 212, 1575, 1259, 840
Yazu, 4834, 840, 315, 179
Suzuki, 618, 1575, 1146, 525
Akou, 21, 4725, 3129, 1575
Kochi, 28, 2100, 874, 210
Okoze, 28, 5250, 2635, 210
Ainame, 29, 3150, 621, 315
Hage, 1205, 1890, 383, 210
Konoshiro, 21, 2100, 1100, 105
Sayori, 12, 5250, 3916, 1260
Anago, 1637, 2625, 1721, 630
Mebaru, 330, 4200, 1680, 315
Tachiuo, 1211, 2310, 930, 105
Hamo, 1622, 3938, 592, 53
Karei, 41, 6300, 3178, 525
Hirame, 540, 2100, 1496, 210
Aji, 5149, 3150, 789, 210
Saba, 5413, 3675, 625, 53
Ika, 2414, 2888, 1083, 105
Tako, 1844, 1890, 1180, 525
Ebi, 560, 7350, 1420, 525
KurumaEbi, 222, 8925, 6508, 2625
Koiwashi, 1316, 1155, 617, 328
ObaIwashi, 2716, 499, 267, 175
Mentai, 678, 1155, 859, 394
Hamachi, 4772, 998, 632, 105
Isaki, 268, 2625, 1465, 840

最初にするべきことは、やはり相関係数を出すこと。
話を解りやすく進めるため、今回は「中値」と「安値」で見ていくことにする。
この二つの変数をそれぞれ、x1x2という変数に入れ直して、相関係数を見る。

x1 <- X$Naka
x2 <- X$Yasu
cor(x1, x2)

実行結果は以下のようになる。

> x1 <- X$Naka
> x2 <- X$Yasu
> cor(x1, x2)
[1] 0.8704575

相関係数の符号はプラスで、値は、0.81に近い値になっている。
ところで、相関係数の算出式はどようなものであったか?
今一度、思い出してみる。たしか、以下のような式だった。



要するに、分母が二つの変数の標準偏差に相当し分子が共分散となっている。
二つの変数の値が全く同じ場合は、共分散の部分が分散となり、
結果として、相関係数が「1になるのであった。忘れた人は、もう一度復習を。

さて、この式は、二つ変数の実際のバラツキの連動性を、
理屈上のバラツキの連動性で基準化しているということになる。
見方を変えると、この式は二つの変数を対等に見ているとも言える。

はて?「対等」ということはどのようなことか?
また、逆に「対等でない」ということはどのようなことか?

要するに、一方の変数のバラツキから、
もう一方の変数のバラツキを「説明する」という考え方もできる。



β の上に「^(ハット)」が付いているが、これは推定値という意味
数学の世界では慣例的に、ある数式によって導かれた値をこのように表す。

なお、相関係数の式では、二つの標準偏差を使うことで、
分散の値に揃えようとしていたのだけれど、ここでは一方のみで良いので、
分母の部分は分散に相当するようになっている。分子はそのまま。

さて、この式は、一体何を表しているか?

要するに、二つの共分散を一方の変数のみで基準化すると、
もう一方の「理屈上の分散を推定できる、と考えているのである。

さらに突き詰めて考えてみると、一方の値が分かれば、
理屈上のもう一方の値の「増減率」を推定することができる。

さて、言葉で考えても解りにくいので、
とりあえず、散布図を描きながら考えてみることにする。
# まず、推定値を求める
B <- cov(x1,x2)/var(x1)

# xの値の範囲を最低値と最大値で決める。
x1.range <- seq(min(x1),max(x1))

# x値の範囲に対応する推定値を計算する。
x1.estim <- B*x1.range

# まずは、普通に散布図を描く
plot(x1, x2, ylim=c(0,max(x2)), xaxt="n", yaxt="n", xlab="", ylab="")

# 図を重ね合わすためのオマジナイ。
par(new=TRUE)

# 推定され得る値を直線で表す
plot(x1.range, x1.estim, col="red", type="l", ylim=c(0,max(x2)), xlab="", ylab="")

# x軸とy軸の軸ラベルを付ける
title(xlab="x1", ylab="x2")

# 変数x2の平均(横線)と変数x1の平均(縦線)を描く
abline(h=mean(x2), v=mean(x1), lty=2, col="blue")

このように、相関係数とは異なった方法で二つの変数の関係を表すことができる。
今回の例の場合は、「それとなく」線と散布図がフィットしているように見えるが、
この赤い線によって推定されているのは、二つの変数の増減率の関係を表す「傾きだけ。

値そのものを推定したいのであれば、二つの値の中心に原点を合わす必要がある。
要するに、上の図で言うと、二つの青い破線の交点に赤い線を通すようにすれば良い。
どのようにすれば良いかは、頭で考えても解るハズ。中学校までの知識で十分。

うん?という人のために補足説明。

推定の赤い線が、x1平均値の時に、x2平均値を取るようすれば良いので、
以下のような補正項を考えれば良いことになる。



x1の推定された値が、実際の平均よりもy軸方向にズレているのが問題なので、
そのズレの部分を計算しているだけ。大して難しい話では無い。
さて、この部分を先ほどの増減率の式と一緒にしてみる。



つまり、最終的にこのようになる。推定されるべきx2の値は、
補正項の部分とx1の推定された増減率の値を足したものからなる。

これを「単回帰」と呼ぶ。

実は、このブログでは、相関係数から単回帰についての説明を導いたが、
歴史的には、相関係数よりも単回帰の方が先に発見された。
ちなみに、発見した人物はナイチンゲールの親戚でもあるフランシス・ゴールトン

少し話が逸れたが、この回帰直線相関係数を同時に表示すると、
さらに、二つの変数間の関係がよく解るのである。
そのため、散布図を描いた際に、回帰直線を一緒に描くことがよくある。

そう言えば、先程の図は、まだ補正項の部分が修正されていなかった。
折角なので、補正した図を作成してみる。
# まず、推定値を求める
B <- cov(x1,x2)/var(x1)

# まず、推定値を求める
B0 <- mean(x2) - B*mean(x1)

# xの値の範囲を最低値と最大値で決める。
x1.range <- seq(min(x1),max(x1))

# x値の範囲に対応する推定値を計算する。
x1.estim <- B0 + B*x1.range

# まずは、普通に散布図を描く
plot(x1, x2, ylim=c(0,max(x2)), xaxt="n", yaxt="n", xlab="", ylab="")

# 図を重ね合わすためのオマジナイ。
par(new=TRUE)

# 推定され得る値を直線で表す
plot(x1.range, x1.estim, col="red", type="l", ylim=c(0,max(x2)), xlab="", ylab="")

# x軸とy軸の軸ラベルを付ける
title(xlab="x1", ylab="x2")

# 変数x2の平均(横線)と変数x1の平均(縦線)を描く
abline(h=mean(x2), v=mean(x1), lty=2, col="blue")

今回は、仕組みを理解するために手計算で回帰直線を引いたが、
もちろん、Rには簡単に回帰分析を行うための関数がある。以下は、実行結果。

> lm(x2 ~ x1)
Call:
lm(formula = x2 ~ x1)
Coefficients:
(Intercept)           x1
   -69.6845       0.3637

Rでは、lm()関数というものを用いる。括弧の中身が解りにくいかもしれない。
これは、モデル式と呼ばれるもので、「」の左側に求めたい変数を置き、
右側に求めたい値を「説明するための変数」を置く。

統計学では、説明するための変数のことを「説明変数」と呼び、
その説明変数によって説明される変数のことを「目的変数」と呼ぶ。

ちなみに、

数学では、説明するための変数のことを「独立変数」と呼び、
その説明変数によって説明される変数のことを「従属変数」と呼ぶ。

文系出身者の多くには、二つの分野が同じように見えるかもしれないが、
前者は工学系の分野で、後者は理学系の分野
論理の立て方が異なるので、与えられる用語は(記号も)微妙に異なる。

実際には同じなのだけれど、用語の使い方に惑わされないように。

さてさて。出力結果で見ないといけないのは、Coefficientsと書かれている部分。
(Intercept) というのは、日本語では「切片」と訳される。要するに、 の部分。
そして、x1と書かれている部分が「傾き」の係数部分になる。

lm()関数は、単回帰だけでなく、重回帰分析でも使用する関数なので、
今の間に覚えておいた方が良いかもしれない。


ふむ。確かに、推定するための線があった方が散布図は見易そうである。
ところで、この推定された値は、実際にはどの程度上手く当てはまっているのだろうか?
これを検討するには、推定された値と、実際の値とのズレを分析する必要があるのだが…。

今回も、話が少し長くなったので、この「ズレ」の検討については次回にしよう。

2012/08/24

文系のための「数の可視化」(1)

ナイチンゲールは、統計学に何をもたらしたのか、今一度、思い出してみる。
重要なことは、統計データの可視化であった。

彼女は、戦場の劣悪な衛生状況が原因で、多くの兵士が感染症で死亡し、
また、その死者数は、戦闘による死者数よりも多いことを知った。
彼女は、その状況を可視化して伝え、可視化された統計データを用いて、
戦場における衛生状況の改善の必要性を説き、当時の政府を説得した。
これが功を奏し、彼女らの提案は受け入れられ、戦場の衛生状況は改善された。

さて、可視化の意味とは何であるか?一言で表現するならば、全体像の抽象化

実は、文章による情報伝達は、実際には効率が悪い
文章は、最初の一文字目から最後の一文字目までが、
一次元の線上に並んでいるようなものであるので、
最後まで読まないと、全体像を把握することができない。

5ページ程度なら良いが、1000ページを超えるレポートならどうか?

なるほど、それでも、文章の方が良いという人もいるかもしれない。
しかしながら、現在のように、情報が超高速で流通している時代に、
そして何より、少子化による極度の人材不足の時代に、
そのような悠長なことは言ってられない。

我々が使える時間は、明らかに、20年前よりも減っている。
「下を探せ」とか「次の世代」などと言う人間には、
理解できないことであろう。

まぁ、良い。とにかく、文章による情報伝達に対して、
グラフィカルな情報は、平面的に表すことができる
文章が一次元なら、グラフィックは二次元の情報と言える。

つまり、全体像を最初から、一瞬で示すことができる。
これこそが、可視化のメリットである。

では、具体的にはどのような方法があるのか?
代表的なものを上げるならば、以下の五種類が一般的であろうか。
  • 棒グラフ
  • 円グラフ
  • 折れ線グラフ
  • 散布図
  • ヒストグラム
上記の5種類のグラフは、用途に合わせて細分化されているし、上記以外のグラフもある。
目的と扱っているデータの種類によって、表現方法は無限に存在し得る

さて、色々と前フリが長くなってしまったが、
グラフを描く上で、最低限、知っておくべきことがある。

グラフィカルな表現は、相手を説得する上で協力なツールとなるが、
その一方で、相手に間違った情報を植え付けるということもあり得る。
実際に、そのようにして、世論を操作する方法もある。
そもそも、グラフを描く上で必要なことは何か、改めて考える余地がある。

さて、今回は、Rでのグラフの描き方について整理してみる。今回は基本だけ。
「共分散」の話で登場した散布図を例に実演してみる。

まずは、読み込み。いつもの作業。


# Windows と Linux の人は次のコマンド
X <- read.table("clipboard", sep=",", header=TRUE)

# Mac の人は次のコマンド
X <- read.table(pipe("pbpaste") , sep=",", header=TRUE)

以下が、そのデータ。前回のものを再掲しただけ。
データの説明に関しては、平均の話を参照。

Oroshi_Kg, Taka, Naka, Yasu
Tai, 1336, 3150, 1217, 53
Chinu, 226, 630, 513, 105
Sawara, 212, 1575, 1259, 840
Yazu, 4834, 840, 315, 179
Suzuki, 618, 1575, 1146, 525
Akou, 21, 4725, 3129, 1575
Kochi, 28, 2100, 874, 210
Okoze, 28, 5250, 2635, 210
Ainame, 29, 3150, 621, 315
Hage, 1205, 1890, 383, 210
Konoshiro, 21, 2100, 1100, 105
Sayori, 12, 5250, 3916, 1260
Anago, 1637, 2625, 1721, 630
Mebaru, 330, 4200, 1680, 315
Tachiuo, 1211, 2310, 930, 105
Hamo, 1622, 3938, 592, 53
Karei, 41, 6300, 3178, 525
Hirame, 540, 2100, 1496, 210
Aji, 5149, 3150, 789, 210
Saba, 5413, 3675, 625, 53
Ika, 2414, 2888, 1083, 105
Tako, 1844, 1890, 1180, 525
Ebi, 560, 7350, 1420, 525
KurumaEbi, 222, 8925, 6508, 2625
Koiwashi, 1316, 1155, 617, 328
ObaIwashi, 2716, 499, 267, 175
Mentai, 678, 1155, 859, 394
Hamachi, 4772, 998, 632, 105
Isaki, 268, 2625, 1465, 840

そして、とりあず、散布図を作成する。散布図とは何か?
たしか、共分散あるいは相関係数に対応する可視化方法であり、
異なる二つの変数間の関係を観察するためのものであった。

卸売数量」と「高値」の関係の基本的なプロットは以下の通り。
plot(X$Oroshi_Kg,X$Taka)

前回までは、詳しく述べなかったが、実は、グラフとして色々と不足している。
まず、このグラフだけでは何のグラフか解らない。やはり、タイトルは必要だろう。

plot(X$Oroshi_Kg,X$Taka)
title(main="卸売数量と高値の関係")

グラフにタイトルが入った。このように、title()関数で、
図にタイトルを入れることができる。
なお、タイトルに関しては、plot()関数でも指定できる。以下は、その例。
結果は、全く同じである。一行で書くのであれば、二つ目の方法となる。


plot(X$Oroshi_Kg, X$Taka, main="卸売数量と高値の関係") 

次に、「x軸ラベル」と「y軸ラベル」を観察してみる。
自動的に付いているが、データを指定した際の名前が自動的に付されていて、
できれば、もう少し、気の利いたラベルを付したい。そこで、以下のようにする。

plot(X$Oroshi_Kg,X$Taka, xlab="卸売数量", ylab="高値", main="卸売数量と高値の関係")

今回は、xlab ylab の二つのパラメータを設定した。これで、見やすい図になってきた。
以上が基本となる。図を作成する際には、タイトルと軸ラベルの調整をすることは重要。

さて、この状態でも、それなりに見やすい図になっているが、
より、見た目を美しい図にするには、どうすれば良いか。以下はセンスの問題。
plot(X$Oroshi_Kg, X$Taka, xlab="卸売数量", ylab="高値", xlim=c(0, 6000), las=1, pch=20, col="red")
grid(col="gray")
text(X$Oroshi_Kg+200, X$Taka + 150, labels=rownames(X), cex=0.5)
title(main="卸売数量と高値の関係")

この例では、plot()関数において、次のパラメータを設定している。
  • xlab: x軸ラベルの名前
  • ylab: y軸ラベルの名前
  • xlim: x軸における値の下限と上限(ylimもある)。
  • las: 軸の値ラベルの描画方向(0~4)で指定。
  • pch: 布置するポイントの記号。20番は塗潰しの丸。
  • col: 布置するポイントの色。今回は、赤色を指定。
さらに、この例では、グリッド線を入れるためのgrid()関数や、
対象名をグラフに布置するめのtext() 関数を用いている。

text()関数に関しては、ポイント記号と重ならないように、
実際のデータよりも、x方向に200、y方向に150、ずらしている。

Rでは、グラフの描き方は非常に多彩であり、様々な図を作成することができる。
もしも、本気でRを用いたいという人は、Rグラフィクスの本を持っていても損は無いだろう。