文系のための「データのバラツキ」（２）

ようやく、データ分析の入り口に辿り着いた。
ここまでの話は、基礎的な話が多かった。

今日の話は、データの「要約」について。

元のデータを分析する前に、データの全容を把握することは重要である。
本来、データの可視化は、全体的な状況を把握してから行うべきであって、
可視化してから考えるのではない。

まず、基本的な情報として考えられるものは何だろうか？
これまでに、やってきた情報を整理すると、

• 行数＝対象の数
• 列数＝属性（変数）の総数
• 各変数の平均と中央値
• 各変数の分散と標準偏差
• 各変数間の共分散と相関係数

と言ったところか。

これだけの情報が得られるだけで、様々なことが解るのだが、
今回は、これらの指標に加えて、「分位数（quantile）」について考えてみる。

分位数は、データのバラツキを、分散とは違った方法で観察する方法とも言える。

とにかく、説明するよりも、手を動かした方が理解しやすい。
とりあえず、いつものように...

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

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

データを読み込む。

Aji,Saba,Buri,Hirame,Karei,Kasago,Mebaru,Kawahagi,Tai,Isaki,Ika,Tako
Hokkaido,0,2,18,12,1210,7,452,0,0,0,660,3
Aomori,0,0,52,14,111,2,66,0,77,0,25,7
Iwate,0,3,2,0,14,0,230,0,0,0,12,2
Miyagi,1,75,25,33,2140,0,116,0,0,0,1,1
Akita,21,1,33,4,0,1,76,0,55,0,4,5
Yamagata,110,0,19,6,0,0,14,0,24,0,7,0
Fukushima,0,0,0,0,72,0,60,0,0,0,0,0
Iabaraki,14,57,39,45,36,20,112,1,79,0,109,46
Chiba,381,96,92,100,1,29,65,34,154,135,328,2
Tokyo,38,6,0,0,3,1,0,2,1,0,0,2
Kanagawa,683,185,823,1,1,32,56,29,168,21,256,11
Nigata,194,25,56,5,2,1,96,1,161,0,483,0
Toyama,1,0,21,0,0,0,0,0,3,0,9,0
Ishikawa,88,17,25,6,0,5,36,5,32,0,0,0
Fukui,160,5,307,4,0,6,3,23,170,1,157,0
Shizuoka,56,109,58,11,0,65,28,5,28,205,8,0
Aichi,20,67,98,1,1,20,8,3,150,22,5,1
Mie,19,15,106,0,0,18,13,1,210,77,61,0
Kyoto,32,5,81,1,5,15,6,3,91,4,9,0
Osaka,32,0,2,0,0,5,3,0,14,0,4,0
Hyogo,93,23,81,9,0,125,68,3,20,0,15,49
Wakayama,22,80,99,6,0,6,15,4,116,121,59,0
Tottori,16,0,3,0,0,0,2,0,29,14,41,0
Shimane,24,3,15,7,0,11,6,17,25,16,155,0
Okayama,1,0,0,0,0,34,2,0,40,0,1,11
Hiroshima,64,7,1,0,1,45,48,1,17,0,10,24
Yamaguchi,274,27,366,49,3,52,94,66,308,16,1,0
Tokushima,16,10,44,5,0,21,14,1,29,6,96,0
Kagawa,0,0,0,0,0,0,0,0,2,0,0,0
Ehime,58,0,9,0,0,25,9,5,55,2,0,2
Kochi,0,2,1,0,0,0,0,0,2,7,0,0
Fukuoka,0,0,0,0,0,0,0,0,0,48,82,0
Nagasaki,0,0,384,22,0,131,0,0,306,155,0,0
Kumamoto,55,0,3,0,0,53,0,0,123,2,0,0
Oita,544,97,4,17,0,62,2,217,75,12,32,0
Miyazaki,49,27,3,0,0,0,0,1,35,63,0,0
Kagoshima,109,52,130,1,0,56,0,1,140,23,7,25
Okinawa,20,0,21,0,0,0,0,4,7,3,0,0

このデータは、ヒストグラムの話で用いたデータと同じ。
農林水産省の「都道府県別魚種別遊漁採捕量」のデータを編集して作ったデータ。
架空とまでは言わないが、現実の分析には適さない。

では、早速。分位数を出してみる。まずは、鯵の分位数を計算。

quantile(X$Aji, probs=seq(0,1,by=0.25))

ここでは、quantile()関数を用いて分位数を計算している。
quantile は、日本語では「分位数」を意味する。

この何分割するかを決めているのが、probs パラメータである。
四分位の場合には、probs=seq(0,1,by=0.25) となっているということは、
0〜1までを、0.25刻みで「四等分」している。

このように、4つに分ける分位数のことを「四分位（quartiles）」と呼ぶ。

つまり、0.00, 0.25, 0.50, 0.75, 1.00 というベクトルによって指定しているわけである。
ここで、注意するべきは、0番目の値が存在すること。したがって、
四分位の場合は、5つの値が返ってくる。

中央値の話でも述べたが、様々な分位数が存在し、それぞれに名前が付いている。

• 二分位（median）＝中央値
• 三分位（tertiles） 
• 四分位（quartiles）
• 五分位（quintiles）
• 六分位（sextiles）
• 十分位（deciles） 
• 十二分位（duo-deciles）
• 二十分位（vigintiles）
• 百分位（percentiles）
• 千分位（permiles）

Rの quantile()関数で、これらを計算することは難しくない。

quantile(X$Aji, probs=seq(0,1,by=0.5))    # 中央値
quantile(X$Aji, probs=seq(0,1,by=0.2))    # 五分位
quantile(X$Aji, probs=seq(0,1,by=0.1))    # 十分位

実行結果は以下の通り。

> quantile(X$Aji, probs=seq(0,1,by=0.5))    # 中央値
  0%  50% 100%
   0   23  683
> quantile(X$Aji, probs=seq(0,1,by=0.2))    # 五分位
   0%   20%   40%   60%   80%  100%
  0.0   0.4  19.8  40.2 102.6 683.0
> quantile(X$Aji, probs=seq(0,1,by=0.1))    # 十分位
   0%   10%   20%   30%   40%   50%   60%   70%   80%   90%  100%
  0.0   0.0   0.4  14.2  19.8  23.0  40.2  57.8 102.6 218.0 683.0

色々とあるが、四分位が最も一般的であるので、
quantile()関数の既定は四分位。したがって、パラメータ指定も不要

quantile(X$Aji)

パラメータ指定無しで、このようにすると、以下のように四分位が計算される。

> quantile(X$Aji)
  0%  25%  50%  75% 100%
   0    1   23   82  683

ここで、四分位の考え方を整理してみる。

• 全体の下から0/4（0%）番目の数：第0四分位数（最小値）
• 全体の下から1/4（25%）番目の数：第1四分位数
• 全体の下から2/4（50%）番目の数：第2四分位数（中央値）
• 全体の下から3/4（75%）番目の数：第3四分位数
• 全体の下から4/4（100%）番目の数：第4四分位数（最大値）

となる。ここで、第1四分位数〜第3四分位数の範囲のことを四分位範囲と呼ぶ。

これは、データの下から25%番目目と75%番目の間の範囲である。
別の言い方をするならば、中央値を中心に全データの50%が収まる範囲である。

英語では、Inter-Quantile Range と呼び、IQRと表記する。
Rには、IQR()関数が存在し、四分位範囲も容易に計算できる。

IQR(X$Aji)

実行結果は、以下の通り。

> IQR(X$Aji)
[1] 81

今回のデータの場合は、1〜82までの範囲なので、
四分位範囲は、「82-1=81」となり、正常に計算されている。

やはり、数理的な話も必要であろう。ここからは、数学的な話。

実は、色々と前提となる知識を必要とするので、
色々と正しくない表記も含まれる。表記方法は変則的。厳密な説明ではない。
Wikipedia の解説では、分位数を関数として定義し、実数空間への投影で説明。

意味の解る人は、Wikipedia へ。
意味の解らない人は、以下の説明へ。

まず、以下のようなベクトルがあったとして、



これを昇順に並べ変えたものを次のように定義する。



総数を「n」とし、求めたい分位数を「q」とすると、
次のように表すことができる。



うん？添字の部分が...何番目？落ち着いて考えてみる。

通常は、1から数えるが、コンピュータでは0番目を基準とする場合が多い。
そこで、「(n-1)q」のところで、一番目を「0」からの順番に書き直している。

少し、混乱するかも。要するに、元のデータが以下のようになっていて、



これを、「0」からの順番で書きなおし、



「1」ずらした分を後で、正しい位置に戻すために「+1」としている。
少々、回りくどい表現であるが仕方あるまい。

さて、これで計算できるかと言うと、実は以下のような条件が付随する。
添字部分が、「(n-1)q+1」だと解りにくいので、とりあえず「t」と置く。



見慣れない記号が一杯出てくる。とりあえず、順番に。

まず、「」の記号は、ある集合に含まれていることを表す記号で、
逆が、「」 の記号。つまり、ある集合に含まれていない状況を表す。

では、一体何の集合なのかというと、 という記号が教えてくれている。
これは、Natural の「N」だから、自然数を表している。自然数は1以上の整数。

要するに、上の条件というのは、「整数」だったらそのままで良いけれど、
「小数点」の値が出てきた場合には、再配分し直しましょう、と言っている。

では、早速、この式に当てはめてRで確認する。
仮に、第1四分位を求めたいのであれば、q = 1/4 = 0.25 となる。

# まずは、データをソートする。
x <- sort(X$Aji)

# 第1四分位数を算出する
(t <- (length(x)-1)*0.25+1)

実行結果は以下の通り。

> (t <- (length(x)-1)*0.25+1)
[1] 10.25

なるほど。見事に、自然数では無い。
ということで、厄介な自然数で無い場合の処理。

((ceiling(t)-t)*x[floor(t)]+(t-floor(t))*x[ceiling(t)])

そして、その結果は、以下の通り。

> ((ceiling(t)-t)*x[floor(t)]+(t-floor(t))*x[ceiling(t)])
[1] 1

確かに、合っているようだが、「1」というのは特殊すぎて自身が無い。
そこで、第1四分位〜第4四分位を連続して計算してみる。

x <- sort(X$Aji)    # まずは、データをソートし、ベクトルxに格納
n <- length(x)        # ベクトルxの数を数える

# For文による繰り返し処理
for(q in seq(0, 1, by=0.25)){
    # 四分位数を算出する
    t <- (n-1)*q+1

    # 自然数かをチェックする。
    if (t > 0 && identical(round(t),t)==TRUE) {
        # 自然数の場合はそのまま
        print(x[t])
    } else {
        # 自然数で無い場合は調整する
        print((ceiling(t)-t)*x[floor(t)]+(t-floor(t))*x[ceiling(t)])
    }
}

さて、ここでは、length()関数とif文、identical()関数が初登場である。

まず、length()関数は、ベクトルに含まれるスカラーの個数を返す関数である。
行列のnrow()関数やnrow()関数に似ているが、行列に対してlength()関数を使うと、
その行列に含まれる全成分の総数を返す。

 identical()関数は、二つの値が同じあるかを確認するための関数であり、
同じ場合には、「TRUE」を返し、異なる場合には「FALSE」を返す。
ここでは、四捨五入した値と元の数字を比較するために用いている。

if文は、for文と同様に、プログラミングの基本となる構文であり、
ある処理の条件分岐を定義する。

ここでは、「t」が自然数であるか、否か、の分岐を示している。
括弧内の処理が、TRUEの場合は、そのままの値を出力し、
それ以外の場合（else）の時には、値を調整して出力している。

この分岐処理の内部では、二つの条件が含まれている。
一つ目は、t > 0 という条件（1以上の正の整数であるという条件）
二つ目は、t が小数点を含んでいるか否か？という条件である。

この二つの条件は、「&&」で結ばれている。これをAND演算と呼ぶ。
すなわち、両方の条件が揃った場合を「真」とする条件が与えられている。

さて、この処理の実行結果は以下の通り。

[1] 0
[1] 1
[1] 23
[1] 82
[1] 683

なるほど、Rの quantile()関数 の結果と見事に同じである。
以上が、分位数の算出方法である。

ところで、今回は、分位数の算出方法を、



で説明した。Excel や R の既定（type 7）はこの方法に従っている。

しかしながら、分位数の分割方法は9通りも存在し、
この辺りは実に厄介な話である。実を言うと、良く理解できない。

この四分位に対して、ヒンジという考え方がある。
これは、四分位と異なり、考え方は至ってシンプル。

中央値より上側の中央値を上側ヒンジと呼び、
中央値より下側の中央値を下側ヒンジと呼ぶ。

R では、fivenum()関数によって計算することができる。

fivenum(X$Aji)

実行結果は以下の通り。

> fivenum(X$Aji)
[1]   0   1  23  88 683

この方法は、五数要約と呼ばれ、fivenum()関数は、
最小値、下側ヒンジ、中央値、上側ヒンジ、最大値、の5つの値を返す。

実用上は、四分位と五数要約を意識することは無いが、
箱ひげ図の場合には、五数要約を用いる。

最後に、今回のデータの四分位をまとめて出力してみる。
これは、summary()関数を用いることで、簡単に計算できる。

summary(X)

行列データから、直接的計算してくれる。非常に楽。
よく見ると、四分位だけではなく、平均値も計算されている。

> summary(X)
      Aji              Saba             Buri           Hirame       
 Min.   :  0.00   Min.   :  0.00   Min.   :  0.0   Min.   :  0.000  
 1st Qu.:  1.00   1st Qu.:  0.00   1st Qu.:  3.0   1st Qu.:  0.000  
 Median : 23.00   Median :  5.00   Median : 23.0   Median :  1.000  
 Mean   : 84.08   Mean   : 26.21   Mean   : 79.5   Mean   :  9.447  
 3rd Qu.: 82.00   3rd Qu.: 27.00   3rd Qu.: 81.0   3rd Qu.:  8.500  
 Max.   :683.00   Max.   :185.00   Max.   :823.0   Max.   :100.000  

     Karei             Kasago           Mebaru          Kawahagi     
 Min.   :   0.00   Min.   :  0.00   Min.   :  0.00   Min.   :  0.00  
 1st Qu.:   0.00   1st Qu.:  0.00   1st Qu.:  0.50   1st Qu.:  0.00  
 Median :   0.00   Median :  6.50   Median : 11.00   Median :  1.00  
 Mean   :  94.74   Mean   : 22.32   Mean   : 44.74   Mean   : 11.24  
 3rd Qu.:   1.75   3rd Qu.: 31.25   3rd Qu.: 63.75   3rd Qu.:  4.00  
 Max.   :2140.00   Max.   :131.00   Max.   :452.00   Max.   :217.00  

      Tai             Isaki             Ika              Tako       
 Min.   :  0.00   Min.   :  0.00   Min.   :  0.00   Min.   : 0.000  
 1st Qu.:  8.75   1st Qu.:  0.00   1st Qu.:  0.25   1st Qu.: 0.000  
 Median : 33.50   Median :  2.00   Median :  8.50   Median : 0.000  
 Mean   : 72.26   Mean   : 25.08   Mean   : 69.39   Mean   : 5.026  
 3rd Qu.:121.25   3rd Qu.: 19.75   3rd Qu.: 60.50   3rd Qu.: 2.000  
 Max.   :308.00   Max.   :205.00   Max.   :660.00   Max.   :49.000  

四分位は、分散と同様にデータ全体のバラツキを観察する上で重要であるが、
この状況を、眺めても、いまいち解りにくい。やはり、標準偏差の方が解かりやすいか？

しかしながら、この状況を可視化することで、非常に重要な参考資料となり得る。
その可視化方法こそが、「箱ひげ図」なのである。この話については、次の話で。

文系のための「データのバラツキ」（１）

今回、少々、長い。画像やデータが多いということもある。
しかし、今回の話は、最重要。理解できないと、先に進めない。

早速、本題に入る。データの「真ん中」には、様々な定義があった。

一般的なものは、平均、中央値、最頻値、など。

「真ん中」とは何か、を突き詰めることは面白い。

しかし、本ブログの目的は、多次元データ解析。先に進む。

さて、そもそも、データの「真ん中」の話をした理由は、

各々の対象が「真ん中」から離れている度合いを知りたいからであった。

なぜ、そのようなことを知りたいのか？

「離れている度合い」によって、対象の「個性の強さを知る」ためである。

ふむ。解ったような、解らないような。では、簡単な質問をする。

対象全体としての「個性」の「バラツキ」は、大きい方が良いか？

それとも、対象全体としての「バラツキ」は、小さい方が良いか？

全体としての「バラツキ」が、小さいということは、

全体として「個性」が「無く」、観察できる視点は限定される。

したがって、答えは、「バラツキ」が大きい方。

何故、そのようなことが言えるのか、

冷静に考えてみれば、簡単に理解することができる。

例えば、教室Aに出席している100人の学生に「将来なりたい職業は？」と質問し、

100人全員が「地方公務員」と回答したとする。

そして、別の教室Bに出席している別の100人に同じ質問をして、

10人が「大学教授」、20人が「弁護士」、20人が「医者」、

30人が「公認会計士」、25人が「起業家」、残り5人が「宇宙飛行士」、

と答えたとする。個性的なクラスである。熱意と野望と夢に満ちている。

できれば、このような個性的なクラスを受け持ちたい。まぁ、良い。次に進む。

さて、教室Aは、全ての学生が「地方公務員」と答えた。夢の無いクラス。
一方、教室Bは、回答が大きく別れた。個性に溢れた夢のあるクラス。

分析対象として、興味が湧いてくるのは明らかに教室Bである。なぜなら、

• 性別の差は存在するのか？
• 年齢が微妙に異なっているのではないか？
• 出身地が影響のではないか？
• 両親の職業が関係している可能性は？

などなど。バラツキが大きい、すなわち、個性が豊かであるということは

新しい問題を発見する上で重要な手がかりとなる。

それゆえに、データ分析では「バラツキ」が重要となる。

かつて、「ゆとり教育」というものがあった。

あの政策が出てきた理由は「個性豊かたな人間を育む」とか、

たしか、そのような理由であった。評価はともかく、「理念」は間違っていない。

とにかく、分析において「バラツキ」というのは、「個性」の象徴であるので、

「バラツキ」の大きいことは良いことなのである。分析者としては。

では、いよいよ、バラツキについて本格的に考えていく。

とりあえず、平均の話で用いた「魚の卸売価格」のデータを用いて考える。

何度もやってきたので、データの読み込み方は大丈夫であろう。

# 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

では、早速、このデータの「個性」を見てみる。

これまでの説明から、最初に、データの「真ん中」を算出し、

各対象の値からの「離れ具合」を求める。とりあえず、以下を順次実行していく。
今回は、とりあえず「平均」を「真ん中」にして考えることにする。
（実際には、中央値の方が適しているが、その議論はいずれ。）

# ステップ1：データの総数（n）を求める。

n <- nrow(X)

ステップ1は、大丈夫だろう。行数を変数nに代入しているだけ。

# ステップ2：卸売数量の平均を求める。転置が必要。

m <- t(colMeans(X))

colMeans()関数は、行列の列平均を算出してくれた。
便利であるが、数値の部分のみを取り出すには「転置」が必要。
これは、数学的な問題ではなく、Rというコンピュータ言語の仕様の問題。

# ステップ3：行列の引き算ができるように、平均値を総数個複製する。

(M <- matrix(rep(m, each=n), nrow=n, ncol=4))

ここでは、新しくrep()関数が登場している。Replicate（複製する）の略。
あるスカラーの値を特定のルールに従って複製することができる。
この場合は、「each = n」つまり「それぞれの値をn回複製する」という意味。
では、変数nには何が入っていたか？確か、行数（対象数）であった。
ここでは、rep()関数は、matrix()関数の中に入っていて、行数はnとなっている。
したがって、元の行列と同じサイズの、平均値が並んだ行列が作られた。

# ステップ4：元の行列から、平均値が複製された行列Mを引き算する。
(A <- X-M)

最後のステップは簡単。単なる行列の引き算。自信が無い人は復習コース。
少し複雑だったかもしれない。慣れれば、大したことは無い。

さて、最後の結果は、新しく作られた変数Aに格納されている。
この行列Aの各成分は、各対象について、平均が引かれたものになっている。
この指標が、各対象の「個性」の強さを表しているものであり、「偏差」と呼ぶ。

ここで少し数式も入れる。既にやったことを数式で表現しているにすぎない。



ここで、対象の数は 1 < i < n 個あり、変数は 1 < j < p 個ある。いつもの「np行列」。

まずは、記号の読み方。「X」の上に「〜」が書いてある。
読み方は「エックス・チルダ」。各成分について「列平均」を引いた行列を表す。
当然、列の数（p個の列）だけ「平均」あるので、列平均を「」のように表現する。

この行列のことを「偏差行列」と呼ぶ。偏差が並んでいる行列。頻繁に登場する。
元の多次元データから、各対象の「個性」だけの行列に作り替えたものである。

では、偏差行列を用いて何ができるのか？試してみる。

この状況では、少々、観察しにくいので、まずは、可視化してみる。可視化は重要。
今回は、結果を美しく見せるためにちょっとした工夫がしてある。説明は後で。

# 卸売数量の可視化

plot(A$Oroshi_Kg, type = "h", xaxt="n", xlab="", main="Oroshi")

abline(h=0)
axis(side=1, at = (1:n), labels = rownames(A), las=2)

# 高値の可視化
plot(A$Taka, type = "h", xaxt="n", xlab="", main="Takane")
abline(h=0)
axis(side=1, at = (1:n), labels = rownames(A), las=2)

# 中値の可視化

plot(A$Naka, type = "h", xaxt="n", xlab="", main="Nakane")
abline(h=0)
axis(side=1, at = (1:n), labels = rownames(A), las=2)

# 安値の可視化

plot(A$Yasu, type = "h", xaxt="n", xlab="", main="Yasune")
abline(h=0)
axis(side=1, at = (1:n), labels = rownames(A), las=2)

まずは、上記のコマンドの説明から。

基本的に、4つの図は同じ方法で作成されている。

plot()関数を使って、基本的な図を書き、

次に、abline()関数を使ってゼロを通る水平線を書いて、

最後に、axis()関数を使って軸ラベルを書いている。

各関数のパラメータを見てみる。

• plot()関数
• type: プロットの種類。棒で表示するので「h」。「histogram」の「h」
• xaxt: この段階では、x軸を描画しないので「n」。「No」 の「n」
• xlab: 今回は、x軸の軸ラベルを表示しないので空白。
• main: 図の上に表示させるグラフタイトル。
• abline()関数
• h: 水平線を通る線を指定。通る場所は「0」。
• axis()関数
• side: ラベルを配置する場所。下に表示するので「1」を指定。
• at: 対象のラベルを配置する位置。等間隔に「1〜総数個」。
• labels: 実際に配置する文字列。行列Aの列名を使うので「rownames(A)」。
• las: ラベルの向き。軸に対して垂直に置くので「2」。

Rでは、様々なパラメータによって図を整える。

上記以外にも、パラメータはあるが、今回は必要なもののみ。

余裕がある人は、Rのヘルプを参照すると良い。詳しい説明がある。

例えば、axis()関数のヘルプを参照するには、

コンソール上で、「?axis()」のように、前にハテナマークを付けるだけ。

次に、結果から何を読み取れるかを考えてみる。

基本的なこととして、平均は丁度0を通る水平線となっていて、

この線よりも上に出ているのは、平均値よりも高く、

この線よりも下に出ているのは、平均値よりも低いことを示している。

とりあえず、「車海老」を観察してみる。

車海老は、卸売数量が平均を下回り、高値、中値、安値が平均よりも高い。

したがって、希少価値が高くて、全体的に価格も高いと言える。

「鯵」は、卸売数量が平均を上回り、高値、仲値、安値ともに平均以下。

つまり、希少価値は低く、全体的に安いと言える。車海老の逆の傾向

面白いのは「チヌ」。卸売数量は、平均よりも低いが、

高値、中値、安値ともに、平均よりも低い。

マダイなど、他のタイ科よりもクセがあるから、

ニーズが限られるのかもしれない。

このように、平均からの離れ具合を見てやると、

それぞれの対象の個性が見えてくる。この感覚が重要である。

今回のデータから、卸売数量と価格の関係を観察したくなるが、
これは、もう少し先の話。もう少し、基礎的な話で我慢する。

さて、今回は、話が長くなる。もう少しの辛抱。
基本的な話は、時して退屈であり、また、そのような話を長々とするのも辛い。

話が逸れた。今度は、全体としての「個性」というのを考えてみる。
全体としての個性というのは、「偏差」を合わせたものと考えることができる。
すでに、「偏差」は変数Aの行列に格納されていた。
これを足しあわせてみるとどうなるか？実行する前に気づいただろうか？
Rでは、列和を次の方法で計算することができる。

colSums(A)

以下は、実行結果。

> colSums(A)
    Oroshi_Kg          Taka          Naka          Yasu
-1.818989e-12  5.002221e-12 -2.501110e-12  0.000000e+00

あれ〜？何かおかしい。若干の誤差があるが、全てが「0」となっている。
なぜか？この場合は、「平均」をデータの「真ん中」にしている。
したがって、偏差の負の側と偏差の正の側を足し合わせると、
結果として、両方を「平（たいら）に均（な）して」しまい、「0」となる。
この方法は、どうやら良くない。では、どうするべきか？

簡単な方法は、最初から、マイナスにならないようにする。
そのような方法は、これまでにやったことがあるだろうか？
思い出してみる。確か、二乗するという方法があった。これを試してみよう。

注意！これから行うのは、行列の計算ではない。行列の「成分の二乗」。

colSums(A^2)

この計算の結果は以下の通り。各列の「二乗和」が計算されている。

> colSums(A^2)
Oroshi_Kg      Taka      Naka      Yasu
 78814194 117904883  49075875   8568607

コンピュータの世界では、べき乗を「^（ハット）」で表す。Rでも同様。
さて、この状況には一つ問題がある。全体の個性を表すことはできているが、
これでは、新しい対象が加わる毎に増えていき、他の状況とは比較できない。
例えば、今回のデータの場合、別の日の「卸売価格」との比較ができない。

どうすれば良いのだったか？各対象の数に合致するように、
対象数の単位で基準化するのだった。つまり、対象数の逆数を掛ける。

colSums(A^2/n)

この計算結果は以下の通り。偏差の二乗の平均。

> colSums(A^2/n)
Oroshi_Kg      Taka      Naka      Yasu
2717730.8 4065685.6 1692271.6  295469.2

この状況も数式を使って一般化してみる。
この式の見方は、既に、平均の話でしている。省略。



「σ」の記号は「シグマ」と呼ぶ。「Σ」の小文字。
二乗和であるため、二乗が付いている。これ重要。

この式は、「母分散」と呼ばれる「バラツキ」の指標である。
統計学では、「母集団」と「標本」という考え方がある。
母集団というのは対象の全てであり、標本はその一部。
多くの場合、「母集団」を特定することは困難なので、
標本をとっていることになる。ふむ。そいうものであるのか。

確かに、今回のデータは、色々な意味で「全体の一部」であった。
ところで、全ての魚介類を一つの市場で扱うことは可能なことなのか？

まぁ、それは良いとして、なるほど、「標本」データである。
実は、この場合、母分散の推定値にはならないことが統計学的に証明されている。
実際の標本の分散は、母分散よりも少し大きくなる。

では、どうするべきか？そこで登場するのが「不偏分散」。
何やら、小難しい言葉が出てきた。ここは、我慢が重要。
それほど複雑ではない。以下のように式を変更する。


変更された部分は一箇所のみ。nの逆数を掛けていた部分を、n-1の逆数で掛ける。
そのような調整をしている。実際に計算してみる。

colSums(A^2/(n-1))

以下は、実行結果。

> colSums(A^2/(n-1))
Oroshi_Kg      Taka      Naka      Yasu
2814792.6 4210888.7 1752709.8  306021.7

実を言うと、R には不偏分散を計算するためのvar()関数がある。
では、上の結果とvar()関数の結果とを比較する。以下を実行。

var(X$Oroshi_Kg)
var(X$Taka)
var(X$Naka)
var(X$Yasu)

以下は、実行結果。計算誤差があるが確かに同じ。
このことから、var()関数は、不偏分散を求める関数であることが解る。

> var(X$Oroshi_Kg)
[1] 2814793
> var(X$Taka)
[1] 4210889
> var(X$Naka)
[1] 1752710
> var(X$Yasu)
[1] 306021.7

ところで、「分散」というのは、偏差の二乗の和となっている。
したがって、全体の「バラツキ」を考えると、実は、二乗分だけ多い。
これは、少々、不都合なことがある。どのような不都合か？
要するに、元の値や、平均値と比較することができない、という不都合。

元に戻すにはどうすれば良いのか、というのがここでの問題。

すでに実験したように、偏差を足しても「0」となる。
そこで、分散の平方根（ルート）を取ることにする。念のため。

　であり　

平方根というのは、二乗して元に戻るような値のことで、
「√（ルート）」という記号によって表したのだった。
Rでは、sqrt()関数というのが用意されているので、これを使う。

sqrt(var(X$Oroshi_Kg))
sqrt(var(X$Taka))
sqrt(var(X$Naka))
sqrt(var(X$Yasu))

そして、以下が実行結果。

> sqrt(var(X$Oroshi_Kg))
[1] 1677.734
> sqrt(var(X$Taka))
[1] 2052.045
> sqrt(var(X$Naka))
[1] 1323.899
> sqrt(var(X$Yasu))
[1] 553.1923

これが、対象数で単位を揃えた偏差の大きさとなる。
つまり、ここで計算された結果は「標準化された偏差」と言える。
一般的には、これを「標準偏差」と呼び、全体としての個性の指標としている。

実は、標準偏差に関しても、Rには関数が用意されている。sd()関数である。

sd(X$Oroshi_Kg)
sd(X$Taka)
sd(X$Naka)
sd(X$Yasu)

いつものように、実行結果を下に。確かに。不偏分散の平方根になっている。

> sd(X$Oroshi_Kg)
[1] 1677.734
> sd(X$Taka)
[1] 2052.045
> sd(X$Naka)
[1] 1323.899
> sd(X$Yasu)
[1] 553.1923

さて、標準偏差に関しても数式で表してみる。慣れが重要。
不偏分散の平方根を取るのだから、難しくはない。



今日の話は長かった。重要なことは以下の3つ。

• 偏差：平均からの離れている程度
• 分散：偏差の二乗の和。全体としてのバラツキの指標
• 標準偏差：分散の平方根をとって、元のデータに合わせた指標。

この3つのことを理解できていれば次の話が解るはず。