R R-console R R Rscript R-console GUI 1

(1)

November 2015

多変量解析入門

(2)

R

超入門

■操作のためのインターフェース

インタラクティブなインターフェース

R-console

の画面を呼び出してコマンドを打ち込む

ターミナルから

R

と打ち込むと，上と同様のキャラクター画面が現れる

ので，コマンドを打ち込む

バッチ処理のためのインターフェース

エディタで

R

のソースファイルを作って，

Rscript

で読み込んで実行

一見むつかしそうだが，

R-console

やターミナルで打ち込むコマンドを

並べるだけなので，この方式に慣れる方がよい．保存性や編集の容易さ

などで計り知れないメリットがある．プログラミングの経験のある人な

ら，この方式を最初から取るほうがよい．

また，

GUI

のメニューで操作するよりも．

1

行ずつ処理を書いていくとい

うやり方のほうが，手続きを言語化することができる分，分かりやすい．

これは大きな利点．

(3)

2

組のデータの相関から

X Y 1 11.04 21.03 2 15.76 24.75 3 17.72 31.28 4 9.15 11.16 5 10.10 18.89 6 12.33 24.25 7 4.20 10.57 8 17.04 33.99 9 10.50 21.01 10 8.36 9.68 5 10 15 20 10 15 20 25 30 35 2変量の散布図 x y

(4)

２変量のデータ

説明変数：

x = [x

₁

, x

₂

, . . . , x

_n

]

目的変数：

y = [y

₁

, y

₂

, . . . , y

_n

]

太字はベクトルであることを表す．つまり，これらのデータは

n

次元の

ベクトル．多変量解析においては．説明変数が

x

₁

, x

₂

, . . .

となる．

行列との積を作る場合には，しばしば縦ベクトルとして扱われる．

x =

         

x

₁

x

₂

...

x

_n          

, y =

         

y

₁

y

₂

...

y

_n          

スペースがもったいないので転置を意味する

T

を付けて次のように書く

こともある．

x = [x

₁

, x

₂

, . . . , x

_n

]

T

(5)

2

変量から作られる統計量

x, y

の２系列のデータから，２つの平均，２つの分散（標準偏差），そし

(6)

平均

¯

x =

1 n

(x

1

+ x

2

+

· · · + x

n

)

¯

y =

1 n

(y

1

+ y

2

+

· · · + y

n

)

だれでも知っている定義だが…

(7)

平均とは何か

?

3

個の数値データ

x

₁

, x

₂

, x

₃

がある．数直線上に並べておいて，動かせる

点

µ

を考える．

3

つの数値データそれぞれとの差の二乗和

S

S = (x

₁

− µ)

2

+ (x

₂

− µ)

2

+ (x

₃

− µ)

2

をとる．

S

が最小になるような

µ

の値は？

µ

µ′

(8)

S

を

µ

で微分してゼロと置いてみると，

dS

dµ

= 2(µ

− x

1

) + 2(µ

− x

2

) + 2(µ

− x

3

) = 0

⇒ µ =

1

3 (x

1

+ x

2

+ x

3

)

平均が現れる！

平均とは，データとの差の二乗和を最小にする値

(9)

分散：平均偏差の二乗の平均

s

2_x

=

1 n

(

(x

₁

− ¯x)

2

+ (x

₂

− ¯x)

2

+

· · · + (x

_n

− ¯x)

2)

s

2_y

=

1 n

(

(y

₁

− ¯y)

2

+ (y

₂

− ¯y)

2

+

· · · + (y

_n

− ¯y)

2)

x

_i

− ¯x

を平均偏差あるいは単に偏差という。

(10)

標準偏差のスケール感覚は大事

n

0 5 10 15 20

0.0

(11)

n

0 5 10 15 20

0.0

(12)

n 0 5 10 15 20 0.0 0.5

n = 1, p = 0.4

の二項分布と平均±標準偏差の領域：

p = 0, 1/2, 1

以外だとはみ出るが，領域の幅は

データの幅を超えることはない．

(13)

共分散

s

_xy

=

1 n

((x

1

− ¯x)(y

1

− ¯y) + · · · + (x

n

− ¯x)(y

n

− ¯y))

(14)

共分散は何を意味するのか？

X

Y

x

−

y

−

I

II

III

IV

I, III

は負の寄与，

II, IV

は正の寄与

(15)

r_xy = −0.769

X

Y

I

II

III

IV

負の寄与のみ

r_xy = 0.737

X

Y

I

II

III

IV

正の寄与のみ

データが右肩上がりに分布しているとき→ 共分散は正

データが右肩上がりに分布しているとき→ 共分散は負

(16)

相関係数は共分散を規格化したもの

(

規格化：最大値を

1

にそろえること

)

r

_xy

=

s

xy

s

_x

s

_y

r

_xy

= 1

完全な正の相関

y = ax + b (a > 0)

全ての点が直線上に乗る

r

_xy

=

−1

完全な負の相関

y = ax + b (a < 0)

r

_xy

≈ 0

無相関

(17)

多変量データはベクトルで表す

X Y 1 11.04 21.03 2 15.76 24.75 3 17.72 31.28 4 9.15 11.16 5 10.10 18.89 6 12.33 24.25 7 4.20 10.57 8 17.04 33.99 9 10.50 21.01 10 8.36 9.68 5 10 15 20 10 15 20 25 30 35 2変量の散布図 x y

散布図では

2

次元の平面上にデータが散らばっている．

しかしこれを

10

次元のベクトルが

2

本ある

というふうに考える．

また，

ベクトルは太字

で

x, y

のように表す．

(18)

物理的なベクトルとは違う

物理的なベクトルは「向き」をもつ量（力，速度，…）

数学的なベクトルは，複数の数の集まり．多変量解析ではデータ数が

次元．

説明変数：

x = [x

₁

, x

₂

, . . . , x

₁₀

]

目的変数：

y = [y

₁

, y

₂

, . . . , y

₁₀

]

しかし，ベクトルの空間的なイメージは理解に役立つ．

要注意！

多変量解析では説明変数が複数になるので，

x

₁

, x

₂

, . . .

と書いて，それぞ

れがベクトルとみなしたい．変量の数が

p

あれば，いわば

x

₁

= [x

_1,1

, x

_1,2

, . . . , x

_10,p

]

こんな感じになる．これは記法上ややこしい．そこで・・・

(19)

平均偏差のほうが扱いやすいので，そちらを使う

生の

(Raw)

データを

x

_R1

, x

_R2

, . . .

のように書いて，

x

を

x = x

_R

− ¯x = [x

_R1

− ¯x, x

_R2

− ¯x, . . .]

のように定義し直す．つまり，特に指定がなくても平均偏差ということ

にする．

y

についても同様．

こうしておくと，分散や共分散がシンプルに書ける．

s

2_x

=

1 n

(x

2 1

+ x

2 2

+

· · · + x

2 n

)

s

_xy

=

1 n

(x

1

y

1

+ x

2

y

2

+

· · · + x

n

y

n

)

(20)

ベクトルのノルム

||x||

2

_{= (x, x) = x}

1

× x

1

+ x

2

× x

2

+

· · · + x

n

× x

n

s

2_x

=

1 n

||x||

2

, s

_x

=

v u u t

1 n

||x||

||x|| : x

のノルム（長さ）

2 つのベクトルの内積

(x, y) = x

₁

× y

₁

+ x

₂

× y

₂

+

· · · + x

_n

× y

_n

s

_xy

=

1 n

(x, y)

(21)

内積の幾何学的な意味

高校数学で出てくるベクトル

⃗a, ⃗b

の内積の定義

(⃗a,⃗b) = ⃗a

の長さ

× ⃗b

の長さ

× cos θ

θ a b O (x1,y1) (x2,y2)

座標成分を使うと，

(⃗a,⃗b) = x

₁

x

₂

+ y

₁

y

₂

(22)

相関係数の幾何学的な意味

r

_xy

=

s

xy

s

_x

s

_y

=

(x, y)

||x|| ||y||

= cos θ

相関係数は，

2

組のデータのベクトルの成す角

度を

θ

とした時の

cos θ

の値

r

_xy

< 0

r

_xy

= 0

r

_xy

> 0

(23)

ベクトルと行列

行列でベクトルを操作する

単位行列

I

：ベクトルを不変にたもつ行列

     

x

y

z

     

=

     

1 0 0

0 1 0

0 0 1

           

x

y

z

     

射影行列：一つの面に「影」を落とすような変換

直交行列：ベクトルを回転させる

(x, y) (x’, y’) O

逆行列

: AA

−1

= I

(24)

連立方程式を扱う

次のような連立一次方程式があるとする．

y

₁

= a

₁₁

x

₁

+ a

₁₂

x

₂

+ b

₁

y

₂

= a

₁₁

x

₁

+ a

₁₂

x

₂

+ b

₂

これを次のように書くと便利．ベクトルを縦ベクトルとして書くことに

注意．

  

y

1

y

₂   

=

  

a

11

a

12

a

₂₁

a

₂₂      

x

1

x

₂   

+

  

b

1

b

₂   

そして一気にまとめて次のように書く．ベクトルは何次元でもよい．

y = Ax + b

x

について解く

(25)

x = A

−1

(y

− b)

上のように解くことで，次の単純な一次方程式の解法と同じ形に持ち込

まれている．

y = ax + b

これから

x = a

−1

(y

− b)

(26)

「独立」とはなにか

２つの確率変数

X, Y

についての大切な関係式

E[X + Y ] = E[X] + E[Y ] (

期待値は「足し算」できる

)

V [X

± Y ] = V [X] + V [Y ] (X, Y

が独立のときだけ

)

標準偏差で考えると，この関係はピタゴラスの定理に似ている．

σ

_X+Y2

= σ

_X2

+ σ

_Y2

確率変数が独立とか，独立でないとかどういうこと？

リレーの選手のタイムで考えてみる．

X, Y, Z

３人の短距離走選手がいる．

Y, Z

の２人は完全に同じ実力．

X

+ Y

の組み合わせと

X + Z

の組み合わせの合計タイムの分布はどうな

るか？

(27)

X σX Y σY X + Y σX+Y σXY ρXY 10.08 0.513 11.50 0.449 21.58 0.962 0.230 0.998 X σX Z σZ X + Z σX+Z σXZ ρXZ 10.08 0.513 11.60 0.369 21.68 0.659 0.018 0.095

X

と

Y

では

σ

_X+Y

= σ

_X

+ σ

_Y

X

と

Z

では

σ

_X+Z2

= σ

_X2

+ σ

_Z2

標準偏差をベクトルとみなすと，独立の時には「直交」

(28)

9.5 10.0 10.5 11.0 11.0 11.2 11.4 11.6 11.8 12.0 12.2 X Y 図1 XとYの調子は完全にシンクロしている 9.5 10.0 10.5 11.0 11.0 11.2 11.4 11.6 11.8 12.0 12.2 X Z 図2 XとYの調子はまったく無関係

(29)

最小二乗法（単回帰分析）

—

重回帰分析の

1

次元版

X

Y

y

=

a

+

b x

(x

1

,y

1

)

(x

2

,y

2

)

(x

i

,y

i

)

(x

i

,a

+

bx

i

)

h

1

h

2

h

i

h

_i

の２乗和が最小になるように，微分を使って

a, b

を決めてやる

(ax + b

の形が多いが，ここでは逆にしてある

)

．

b =

s

xy

s

2_x

, a = ¯

y

− b¯x

(30)

5 10 15 20 10 15 20 25 30 35 xR yR

(31)

重回帰分析

説明変数が２個になった場合

モデル式

y = a + b

₁

x

₁

+ b

₂

x

₂

b

₁

, b

₂

を列（縦）ベクトルで書くと，次のように分散と共分散を成分とす

る行列とベクトルの積で表される．

  

b

1

b

₂   

=

  

s

2 x1

s

x1x2

s

_x₁_x₂

s

2_x 2    −1   

s

x1y

s

_x₂_y   

説明変数同士の共分散

s

_x₁_x₂

が，ここでは影響を与えることに注意．

(32)

一般の場合

(

説明変数が

p

個

)

モデル式

y = a + b

₁

x

₁

+ b

₂

x

₂

+

· · · + b

_p

x

_p

次のようにして，分散と共分散を使って係数

b

₁

, b

₂

, . . .

を求めることがで

きる．

_         

b

₁

b

₂

...

b

_p          

=

         

s

_x₁_x₁

s

_x₁_x₂

· · · s

_x₁_x_p

s

_x₂_x₁

s

_x₂_x₂

· · · s

_x₂_x_p

...

. . .

_...

s

_x_p_x₁

s

_x_p_x₂

· · · s

_x_p_x_p           −1

=

         

s

_x₁_y

s

_x₂_y

...

s

_x_p_y          

なお，ここで

x

₁

等はそれぞれ

n

個のデータからなる．計算量は膨大だ

が，もちろんコンピュータの適切なアプリケーションで高速に処理で

きる．

(33)

R

による単回帰分析

最初の二次元データについて分析を試みる．

データファイル：

xy10.dat

データ：

2

変量

標本の大きさ：

10

X Y 11.04 21.03 15.76 24.75 17.72 31.28 9.15 11.16 10.1 18.89 12.33 24.25 4.2 10.57 17.04 33.99 10.5 21.01 8.36 9.68

DT <- read.table("xy10.dat",header=TRUE)

postscript("Images/lmXY10.ps",

horizontal=FALSE,

height=5.2,width=6,onefile=TRUE)

result = lm(Y ~ X, data = DT)

#

線形回帰を実行

summary(result)

#

結果の数値を出力

plot(Y ~ X, data = DT)

#

散布図出力

(34)

プログラムの意味

DT <- read.table("xy10.dat",header=TRUE)

xy10.dat

というファイルをデータフレームの形式で読み込んで，

DT

と

いう名前のオブジェクトに読み込む．ファイルには

X Y

というヘッダが

あり，

1

列目のデータには

X, 2

列目のデータには

Y

という「名前」が付

与される．

postscript("lmXY10.ps",

horizontal=FALSE,

height=5.2,width=6,onefile=TRUE)

以降の画像出力を

Postscript

形式で書き込む．サイズの単位はインチ

(2.54 cm)

．

result = lm(Y ~ X, data = DT)

(35)

うオブジェクトに代入．

lm

の引数，

Y ~ X, data = DT

は，

DT

オブジェクトの

X

を説明変数，

Y

を目的変数とする意味．

summary(result)

result

の内容，つまり線形回帰の結果のデータを出力．

plot(Y ~ X, data = DT)

データフレーム

DT

の

X

に対して

Y

をプロットして散布図を描く．

abline(result)

線形回帰によって得られた直線を引く．

(36)

このように，

R

のプログラムは，統計処理のためのコマンドを並べていく

だけで，いろいろな分析や描画ができるようになっている．コンピュー

タプログラムというよりも，むしろレシピに近い

（一方で本格的な数値

(37)

結果の出力

Call:

lm(formula = Y ~ X, data = DT) Residuals:

Min 1Q Median 3Q Max

-5.014 -2.754 1.221 2.372 3.491 Coefficients:

Estimate Std. Error t value Pr(>|t|) (Intercept) -0.6092 3.4405 -0.177 0.863859

X 1.8305 0.2800 6.538 0.000181 ***

---Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 3.54 on 8 degrees of freedom Multiple R-squared: 0.8424,Adjusted R-squared: 0.8227 F-statistic: 42.75 on 1 and 8 DF, p-value: 0.000180

(38)

結果を読む

Residuals:

残差の五数要約

Coeﬃcients:

回帰直線の情報：

Intercept

切片

, X

偏回帰係数（傾き）

Estimate

予測値，

Std. Error

標準誤差，

t value t

値，

Pr(>|t|), p

値

(

切片，予測値，

p

値が重要

)

Residual standard error:

残差標準誤差（自由度）

Multiple R-squared:

(

重

)

相関係数の

2

乗

Adjusted R-squared:

(

重

)

相関係数の

2

乗

(

自由度調整済

)

(39)

R

による重回帰分析

データ：ファイル名

BodyScore.txt

Bust West Hip Weight

84 58 87 47 84 59 89 54 86 59 90 50 87 63 94 55 83 60 88 51 83 60 88 50 84 60 90 54 82 60 86 50 82 60 88 52 85 63 90 53 (全30件：ちょっと古い本ではこういうデータを好んで取り上げる傾向があります)

説明変数：

Bust, West, Hip

(40)

データ間の相関行列をチェック

BS <- read.table("BodyScore.txt",header=T)

cor(BS)

##

データ間の相関係数を計算して表示

Bust

West

Hip

Weight

Bust

1.0000000 0.3000223 0.6240064 0.4580098

West

0.3000223 1.0000000 0.5753726 0.6212204

Hip

0.6240064 0.5753726 1.0000000 0.6888909

Weight 0.4580098 0.6212204 0.6888909 1.0000000

一般に，説明変数同士に強い相関があることは好ましくない．ここでは

たかだか

0.6

程度なので問題なし．

(41)

グラフィカルに確認することも有用

Bust 56 57 58 59 60 61 62 63 50 55 60 78 80 82 84 86 88 56 57 58 59 60 61 62 63 West Hip 84 86 88 90 92 94 78 80 82 84 86 88 50 55 60 84 86 88 90 92 94 Weight Bust 56 57 58 59 60 61 62 63 50 55 60 78 80 82 84 86 88 56 57 58 59 60 61 62 63 0.30 West 0.62 0.58 Hip 84 86 88 90 92 94 78 80 82 84 86 88 50 55 60 0.46 _0.62 84 86 88 90 92 94 0.69 Weight

前例の

cor(BS)

の代わりに，

pairs(BS)

とすると，左のように散布図

が表示されて，全体の傾向をつかめる．右はさらに細かい情報を表示さ

せている

(

文献

[2]

第

7

章参照

)

．

(42)

結果を出力する

BS.fit <- lm(Weight ~ Bust + West + Hip, data = BS)

#

説明変数を

+

でつなぐ

summary(BS.fit)

結果：

Residuals:

-4.5320 -1.1779 -0.3508 1.4179 6.4489 Coefficients:

Estimate Std. Error t value Pr(>|t|) (Intercept) -51.8052 19.2769 -2.687 0.0124 * Bust 0.1165 0.2475 0.471 0.6418 West 0.6130 0.2873 2.133 0.0425 * Hip 0.6383 0.2835 2.251 0.0330 * ---Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 2.385 on 26 degrees of freedom

Multiple R-squared: 0.554,Adjusted R-squared: 0.5025 F-statistic: 10.76 on 3 and 26 DF, p-value: 8.859e-05

(43)

主成分分析

多変量データを，なるべく少ない合成関数で要約する．

もしも，入学後の成績に対して数学の点数と物理の点数が同じような寄

与をするとしたら，適当に混ぜ合わせた「数理」の点数を使えばよいの

ではないか？

(44)

p

個の変量ベクトル

x

₁

, x

₂

, . . . , x

_p

がある．これらの重ね合わせとして次

のような

m

個の変量

z

₁

, z

₂

, . . . , z

_m

(m

≤ p)

を構成する．

z

₁

= c

₁₁

x

₁

+ c

₁₂

x

₂

+

· · · + c

_1p

x

_p

z

₂

= c

₂₁

x

₁

+ c

₂₂

x

₂

+

· · · + c

_2p

x

_p

. . .

z

_m

= c

_m1

x

₁

+ c

_m2

x

₂

+

· · · + c

_mp

x

_p

m

の数を少なくでき，かつ

z

_i

(1

≤ i ≤ m)

同士の相関がゼロであるよう

にできれば，これらのベクトルは，

p

よりも小さな次元の空間の中で直

交している．つまり，より少ない情報で記述できることになる．

(45)

数学的な原理

ベクトル

x

₁

, x

₂

, . . . , x

_p

の間の相関係数の行列

R =

         

r

_x₁_x₁

r

_x₁_x₂

· · · r

_x₁_x_p

r

_x₂_x₁

r

_x₂_x₂

· · · r

_x₂_x_p

...

. . .

_...

r

_x_p_x₁

r

_x_p_x₂

· · · r

_x_p_x_p          

の固有値と固有ベクトルを求める．

大きな固有値に属するベクトルほど大きな分散をもつ，主要な成分であ

る．主要な成分を少数のベクトルで構成することができれば，主成分が

うまく掴みだされたことになる．

(46)

例：米国の大気汚染に関わる要因を分析する

SO2 Neg.Temp Manuf Pop Wind Precip Days

Phoenix 10 -70.3 213 582 6 7.05 36 Little Rock 13 -61 91 132 8.2 48.52 100 San Francisco 12 -56.7 453 716 8.7 20.66 67 Denver 17 -51.9 454 515 9 12.95 86 Hartford 56 -49.1 412 158 9 43.37 127 Wilmington 36 -54 80 80 9 40.25 114 Washington 29 -57.3 434 757 9.3 38.89 111 Jacksonville 14 -68.4 136 529 8.8 54.47 116 Miami 10 -75.5 207 335 9 59.8 128 Atlanta 24 -61.5 368 497 9.1 48.34 115 ... ... ... ... ... ... ... ... Charleston 31 -55.2 35 71 6.5 40.75 148 Milwaukee 16 -45.7 569 717 11.8 29.07 123

SO2: 大気中 SO2 の含有量（目的変数），Neg.Temp: 年平均気温（華氏）を負にした値（以下，説明変数），Manuf:

(47)

説明変数同士の散布図や相関を調べてみる

Neg.Temp 0 1000 2500 6 7 8 9 11 40 80 120 160 -75 -65 -55 -45 0 1000 2000 3000 0.19 Manuf 0.063

0.96

Pop 0 1000 2000 3000 6 7 8 9 10 12 0.35 0.24 0.21 Wind 0.39 0.032 0.026 0.013 Precip 10 20 30 40 50 60 -75 -65 -55 -45 40 80 120 160 0.43 0.13 0 1000 2500 0.042 _0.16 10 20 30 40 50 60 0.50 Days 外れ値に注意！

(48)

主成分分析を実行

D <- read.table("usair.txt",header=T)

VD <- D[,-1]

#

最初の変数

SO2

以外の変数

cor(VD)

#

相関係数行列を出力

VD.pc <- princomp(VD,cor=T)

# cor=T

は相関係数行列を使うという意味

summary(VD.pc,loading=T)

#

結果のまとめ．

相関係数行列

Neg.Temp Manuf Pop Wind Precip Days

Neg.Temp 1.00000 0.19004 0.06267 0.34973 -0.38625 0.43024 Manuf 0.19004 1.00000 0.95526 0.23794 -0.03241 0.13182 Pop 0.06267 0.95526 1.00000 0.21264 -0.02611 0.04208 Wind 0.34973 0.23794 0.21264 1.00000 -0.01299 0.16410 Precip -0.38625 -0.03241 -0.02611 -0.01299 1.00000 0.49609 Days 0.43024 0.13182 0.04208 0.16410 0.49609 1.00000

(49)

主成分分析の結果

Importance of components:

Comp.1 Comp.2 Comp.3 Comp.4 Comp.5

Standard deviation 1.4819456 1.2247218 1.1809526 0.8719099 0.33848287 Proportion of Variance 0.3660271 0.2499906 0.2324415 0.1267045 0.01909511 Cumulative Proportion 0.3660271 0.6160177 0.8484592 0.9751637 0.99425879 Comp.6 Standard deviation 0.185599752 Proportion of Variance 0.005741211 Cumulative Proportion 1.000000000

6

つの成分

(components)

の標準偏差，分散の比率（合計で

1

になる），そ

の累計が表示されている．

最初の

3

つの成分だけで全体の

85%

になって

いるので，

3

つの主成分で大体説明できる．

(50)

負荷量

(

元の変量の寄与の割合

)

を見る

Loadings:

Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Neg.Temp -0.330 0.128 0.672 0.306 0.558 0.136 Manuf -0.612 -0.168 -0.273 0.137 0.102 -0.703 Pop -0.578 -0.222 -0.350 0.695 Wind -0.354 0.131 0.297 -0.869 -0.113 Precip 0.623 -0.505 -0.171 0.568 Days -0.238 0.708 0.311 -0.580

•

第

1

主成分：工場の数，人口という社会的要因

•

第

2

主成分：降水量，雨天日数という降雨に関する要因

•

第

3

主成分：負の気温と降水量が逆の寄与

これらのうち，どれが大気汚染に強く関わるだろうか？

(51)

3 つの主成分を散布図で見る

6 つの説明変数を第

1 ，第

2 ，第

3 の主成分

(PC1, PC2, PC3)

に縮約したところで，これ

らの間の関係を散布図で見てみることは有用だ．次のソースを追加することで，必要な

図が得られる．

postscript("PC1vsPC2.ps",horizontal=F, onefile=TRUE)

par(pty = "s")

#

文字の重なりを設定

plot(VD.pc$scores[,1],VD.pc$scores[,2],

#1

列目と

2 列目による散布図

ylim = range(VD.pc$scores[,1]),

#

y

軸の幅を

PC1

に合わせる

xlab = "PC1", ylab = "PC2",

type = "n", lwd = 2)

# type="n"

で枠のみの出力

text(VD.pc$scores[,1], VD.pc$scores[,2],

#

散布図の点をテキストで表示

labels = abbreviate(row.names(D)),cex = 0.7,lwd=2)

dev.off()

##

上のソースは

PC1

と

PC2

の散布図．他の組み合わせは少し修

##

正して繰り返して作る．

(52)

PC1 vs. PC2

の散布図：

Chcg

が外れ値になっている

-6 -4 -2 0 2 -6 -4 -2 0 2 PC1 PC2 Phnx LttR SnFr Dnvr Hrtf Wlmn Wshn Jcks Miam Atln Chcg Indn DsMn Wcht Lsvl NwOr Bltm Dtrt M-SP KnsC St.LOmah Albr Albn Bffl Cncn Clvl Clmb Phld PttsPrvd MmphNshv Dlls Hstn SlLC Nrfl Rchm Sttl Chrl Mlwk

(53)

PC1 vs. PC3

の散布図：

Chcg

が外れ値になっている

-6 -4 -2 0 2 -6 -4 -2 0 2 PC1 PC3 Phnx LttR SnFr Dnvr Hrtf Wlmn Wshn Jcks Miam Atln Chcg Indn DsMn Wcht Lsvl NwOr Bltm Dtrt M-SP KnsC St.L

OmahAlbn Albr Bffl Cncn Clvl Clmb Phld PttsPrvd MmphNshv Dlls Hstn SlLC Nrfl Rchm Sttl Chrl Mlwk

(54)

PC2 vs. PC3

の散布図：

Phnx

が外れ値になっている

-4 -3 -2 -1 0 1 2 -4 -3 -2 -1 0 1 2 PC2 PC3 Phnx LttR SnFr Dnvr Hrtf Wlmn Wshn Jcks Miam Atln Chcg Indn DsMn Wcht Lsvl NwOr Bltm Dtrt M-SP KnsC St.L Omah Albr Albn Bffl Cncn Clvl Clmb Phld Ptts Prvd Mmph Nshv Dlls Hstn SlLC Nrfl Rchm Sttl Chrl Mlwk

(55)

VD.pc

からさらに情報を抽出する

以下を順々にやって，様子を探る．

str(VD.pc) # VD.pc というデータの構造を調べてみる

VD.pc$scores # 6つの成分のスコア

VD.pc$scores[,1:3] # 3つの主成分のスコアだけ表示してみる

Comp.1 Comp.2 Comp.3 Phoenix 2.440096802 -4.19114925 -0.94155229 Little Rock 1.611599761 0.34248684 -0.83970812 San Francisco 0.502073845 -2.25528717 0.22663991 Denver 0.207434109 -1.96320936 1.26621359 Hartford 0.219106349 0.97630584 0.59461329 Wilmington 0.996140738 0.50074082 0.43332129 Washington 0.022928417 -0.05456742 -0.35387289 Jacksonville 1.227849872 0.84912801 -1.87611109 Miami 1.533160553 1.40469861 -2.60660585 Atlanta 0.598994755 0.58723563 -0.99541128 ... Charleston 1.429753320 1.21058365 -0.07944350 Milwaukee -1.391024518 0.15761872 1.69127813

(56)

3 つの主成分は汚染度を説明するか？

汚染度

SO2

に対して主成分１，２，３がどのように効いているかを散布図で確かめ

よう．

postscript("3PCvsSO2.ps",horizontal=F,

height=5.2,width=6, onefile=TRUE)

par(mfrow = c(1,3))

# 3

つのグラフを上の画像内に横に並べて表示

plot(VD.pc$scores[,1], VD[,1], xlab = "PC1", ylab="SO2")

plot(VD.pc$scores[,2], VD[,1], xlab = "PC2", ylab="SO2")

plot(VD.pc$scores[,3], VD[,1], xlab = "PC3", ylab="SO2")

## VD.pc

は主成分分析の結果．その１列目，２列目，３列目を

x

軸にとって

## VD[,1]

つまり元データの１列目すなわち

SO2

の値をプロットしている

(57)

得られた散布図

:

外れ値を除外すると，

PC1

と

SO2

の相関が目立つ

-6 -4 -2 0 2 -75 -70 -65 -60 -55 -50 -45 PC1 SO2 -4 -3 -2 -1 0 1 2 -75 -70 -65 -60 -55 -50 -45 PC2 SO2 -2 -1 0 1 -75 -70 -65 -60 -55 -50 -45 PC3 SO2

(58)

主成分を使った重回帰分析を行なう

## PC1,2,3 を説明変数として，重回帰分析を行なう

pclm <- lm(D$SO2 ~ VD.pc$scores[,1] + VD.pc$scores[,2] + VD.pc$scores[,3]) summary(pclm)

結果：

Residuals:

-36.420 -10.981 -3.184 12.087 61.273 Coefficients:

Estimate Std. Error t value Pr(>|t|)

(Intercept) 30.049 2.907 10.336 1.85e-12 *** VD.pc$scores[, 1] -9.942 1.962 -5.068 1.14e-05 *** VD.pc$scores[, 2] 2.240 2.374 0.943 0.352 VD.pc$scores[, 3] -0.375 2.462 -0.152 0.880 ---Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Residual standard error: 18.62 on 37 degrees of freedom

Multiple R-squared: 0.4182,Adjusted R-squared: 0.371 F-statistic: 8.866 on 3 and 37 DF, p-value: 0.0001473

(59)

因子分析

複数科目の試験の実施データがある．これらの成績は，少数の「共通因子」と，固有の

因子の和として表すことができるのではないだろうか？

表1 スピアマンによる6つのテストの間の相関係数を個々の点数データから計算したもの古典フランス語英語数学音程の弁別音楽的才能古典 1.00 0.83 0.78 0.70 0.66 0.63 フランス語 0.83 1.00 0.67 0.67 0.65 0.57 英語 0.78 0.67 1.00 0.64 0.54 0.51 数学 0.70 0.67 0.64 1.00 0.45 0.51 音程の弁別 0.66 0.65 0.54 0.45 1.00 0.40 音楽的才能 0.63 0.57 0.51 0.51 0.40 1.00

(60)

p

個の変量ベクトル

z

1

, z

2

, . . . , z

p

がある．

(p

はたとえば科目数

)

z

1

= a

11

f

1

+ a

12

f

2

+

· · · + a

1r

f

r

+ u

1

v

1

z

2

= a

21

f

1

+ a

22

f

2

+

· · · + a

2f

f

r

+ u

2

v

2

. . .

z

p

= a

p1

f

1

+ a

p2

f

2

+

· · · + a

pf

f

r

+ u

p

v

p

f

1

, f

2

, . . . :

共通因子

v

1

, v

2

, . . . :

独自因子

a

i,j

:

共通因子パターン（共通因子負荷量）

h

2_j

= a

2_j1

+

· · · + a

2_jr

:

試験

j

の共通性

u

j

:

試験

j

の独自因子負荷量

h

2_j

+ u

2_j

= 1

(

共通因子と独自因子の寄与に分解

)

(61)

いくつかの条件を置く

f

1

, f

2

, . . .

は，すべての

v

i

(i = 1, 2, . . . , )

と直交する

異なる

v

i

, v

j

(i

̸= j)

は直交する

f

1

, f

2

, . . . , v

1

, v

2

, . . .

の分散はすべて

1 f

の直交性については

解き方によって，異なる

f

i

, f

j

(i

̸= j)

が直交する場合（直交解）とそうでない場合

(

斜

交解

)

がある．

(62)

最も簡単なケース

x y z 1 28 29 28 2 18 23 18 3 11 22 16 4 21 23 22 5 26 29 26 6 20 23 22 7 16 22 22 8 14 23 24 9 24 29 24 10 22 27 24

10 人に対する

3 科目の試験で，右表のような点数が得られた

とする．これから相関係数行列を作ると，下のようになる．

x y z x 1.000 0.866 0.787 y 0.866 1.000 0.753 z 0.787 0.753 1.000

共通因数は

1 つだけで，それを

f

とすると，解くべき式は次の形になる．

z

x

= af + u

x

v

x

z

y

= bf + u

y

v

y

z

= cf + u

z

v

z

x

, z

y

, z

z

は平均偏差ベクトル．

(63)

与えられた条件を使って解いていくと，相関係数と

a, b, c

の関係が得られる．

r

xy

= ab, r

yz

= bc, r

xz

= ac

変数

3 つ，方程式の数も

3 つだから，これは解ける．

a =

√

r

xy

r

xz

r

yz

= 0.952

同様にして

b = 0.910, c = 0.828

．

a

2

, b

2

, c

2

が求める共通性．

また，独自因子は

u

2_x

= 0.094, u

2_y

= 0.172, y

_z2

= 0.314

．

このようにして，

どの科目の点数も大きな共通因子１つと独自因子から成る

ことが示さ

れた．

(64)

因子はかならずあるのか

次の相関行列が与えられた．

x y z x 1.00 0.83 0.78 y 0.83 1.00 0.67 z 0.78 0.67 1.00

前の問題と同様に解くと，

a = 1.2, b = 0.7, c = 0.5

u

2_x

=

−0.44, u

2_y

= 0.51, u

2_z

= 0.75

負の分散が現れるので，この解はあり得ない（共通因子を抽出できない）．

共通性が

1 を超えたら結果は使えない．

(65)

因子分析は多義的

最初のデータに因子分析を適用したい．

古典フランス語英語数学音程の弁別音楽的才能古典 1.00 0.83 0.78 0.70 0.66 0.63 フランス語 0.83 1.00 0.67 0.67 0.65 0.57 英語 0.78 0.67 1.00 0.64 0.54 0.51 数学 0.70 0.67 0.64 1.00 0.45 0.51 音程の弁別 0.66 0.65 0.54 0.45 1.00 0.40 音楽的才能 0.63 0.57 0.51 0.51 0.40 1.00

z

1

= a

11

f

1

+ a

12

f

2

+

· · · + a

1r

f

r

+ u

1

v

1

z

2

= a

21

f

1

+ a

22

f

2

+

· · · + a

2f

f

r

+ u

2

v

2

. . .

z

p

= a

p1

f

1

+ a

p2

f

2

+

· · · + a

pf

f

r

+ u

p

v

p

この解には多義性があり，さまざまなアルゴリズムが存在する．

最尤因子分析：確率的なモデルを採用．最もよく使われる

(66)

因子の回転

上の式を行列の形で表現する．

Z = AF + U V

ここで任意の直交行列

R

を使って，次のように

F

を「回転」させてもかまわない．

Z = (AR)(R

−1

F ) + U V

したがって，解は無数に存在する．

さらに

R

は逆行列をもつもの（正則行列）であれば，任意に選ぶこともできる．

因子分析の進め方

適当なアルゴリズムを選んで，因子の数を

1 からひとつずつ増やしていく．

3 程度まで

の少ない因子で説明ができれば結果を採用できる．

(67)

各国の平均余命

次の一覧は，

32 の国

(

地域

)

，年齢，性別ごとの平均余命である．これから平均余命を支

配する因子の数を求めたい．

m0 m25 m50 m75 w0 w25 w50 w75 Algeria 63 51 30 13 67 54 34 15 Cameroon 34 29 13 5 38 32 17 6 Madagascar 38 30 17 7 38 34 20 7 Mauritius 59 42 20 6 64 46 25 8 Reunion 56 38 18 7 62 46 25 10 Seychelles 62 44 24 7 69 50 28 14 South Africa(C) 50 39 20 7 55 43 23 8 South Africa(W) 65 44 22 7 72 50 27 9 Tunisia 56 46 24 11 63 54 33 19 . . . United States (66) 67 45 23 8 74 51 28 10 United States (NW66) 61 40 21 10 67 46 25 11 United States (W66) 68 46 23 8 75 52 29 10 United States (67) 67 45 23 8 74 51 28 10 Argentina 65 46 24 9 71 51 28 10 Chile 59 43 23 10 66 49 27 12 Columbia 58 44 24 9 62 47 25 10 Ecuador 57 46 28 9 60 49 28 11

(68)

因子分析を

1

因子で実行してみる

life <- read.table("../Data/chap4lifeexp.txt",header = T) #

読み込み

life.fa1 <- factanal(life, factors = 1, method = "mle") # 1

因子で分析

(69)

1 因子での結果

Call:

factanal(x = life, factors = 1, method = "mle") Uniquenesses: # 独自性 m0 m25 m50 m75 w0 w25 w50 w75 0.238 0.470 0.399 0.696 0.217 0.005 0.117 0.532 Loadings: # 因子負荷量 Factor1 # 第一因子（ここでは1つしかないのでこれだけ） m0 0.873 m25 0.728 m50 0.776 m75 0.552 w0 0.885 w25 0.998 w50 0.940 w75 0.684 Factor1 SS loadings 5.329 Proportion Var 0.666

Test of the hypothesis that 1 factor is sufficient. # 検定結果に注目! The chi square statistic is 163.11 on 20 degrees of freedom.

(70)

因子を増やしていく

life.fa2 <- factanal(life, factors = 2, method = "mle") life.fa2

life.fa3 <- factanal(life, factors = 3, method = "mle") life.fa3

検定結果だけ引用：

# 因子数２

Test of the hypothesis that 2 factors are sufficient.

The chi square statistic is 45.24 on 13 degrees of freedom. The p-value is 1.91e-05

# 因子数3

The chi square statistic is 6.73 on 7 degrees of freedom.

(71)

因子数

3

についての結果

Uniquenesses:

m0 m25 m50 m75 w0 w25 w50 w75

0.005 0.362 0.066 0.288 0.005 0.011 0.020 0.146 Loadings:

Factor1 Factor2 Factor3

m0 0.964 0.122 0.226 m25 0.646 0.169 0.438 m50 0.430 0.354 0.790 m75 0.525 0.656 w0 0.970 0.217 w25 0.764 0.556 0.310 w50 0.536 0.729 0.401 w75 0.156 0.867 0.280

Factor1 Factor2 Factor3

SS loadings 3.375 2.082 1.640

Proportion Var 0.422 0.260 0.205 Cumulative Var 0.422 0.682 0.887

The chi square statistic is 6.73 on 7 degrees of freedom. The p-value is 0.458

(72)

国ごとの因子得点を求める

ソース：

scores <- factanal(life, factors = 3, method = "mle", scores = "regression")$scores

scores

結果：

Factor1 Factor2 Factor3 Algeria -0.258062561 1.90095771 1.91581631 Cameroon -2.782495791 -0.72340014 -1.84772224 Madagascar -2.806428187 -0.81158820 -0.01210318 Mauritius 0.141004934 -0.29028454 -0.85862443 Reunion -0.196352142 0.47429917 -1.55046466 Seychelles 0.367371307 0.82902375 -0.55214085 South Africa(C) -1.028567629 -0.08065792 -0.65421971 South Africa(W) 0.946193522 0.06400408 -0.91995289 Tunisia -0.862493550 3.59177195 -0.36442148 Canada 1.245304248 0.29564122 -0.27342781 ...

これは見づらい．

(73)

２つの因子についてプロットしてみる

plot(scores[,1], scores[,2], type = "n", xlab = "因子1", ylab = "因子2") text(scores[,1],scores[,2],labels=row.names(life), cex = 1.1, lwd=2) -2 -1 0 1 -2 -1 0 1 2 3 因子1 因子2 Algeria Cameroon Madagascar Mauritius Reunion Seychelles

South Africa(C) South Africa(W) Tunisia Canada Costa Rica Dominican Rep El Salvador Greenland Grenada Guatemala Honduras Jamaica Mexico Nicaragua Panama Trinidad(62)_{Trinidad (67)} United States (66) United States (NW66) United States (W66) United States (67) Argentina Chile Columbia Ecuador

(74)

-2 -1 0 1 -1 0 1 2 因子1 因子3 Algeria Cameroon Madagascar Mauritius Reunion Seychelles South Africa(C) South Africa(W) Tunisia Canada Costa Rica Dominican Rep El Salvador Greenland Grenada Guatemala Honduras Jamaica Mexico Nicaragua Panama Trinidad(62) Trinidad (67) United States (66) United States (NW66)United States (67)United States (W66)

Argentina Chile Columbia Ecuador -2 -1 0 1 2 3 -1 0 1 2 因子2 因子3 Algeria Cameroon Madagascar Mauritius Reunion Seychelles South Africa(C) South Africa(W) Tunisia Canada Costa Rica Dominican Rep El Salvador Greenland Grenada Guatemala Honduras Jamaica Mexico Nicaragua Panama Trinidad(62) Trinidad (67) United States (66) United States (NW66)United States (67)United States (W66)

Argentina Chile Columbia

(75)

因子分析と主成分分析の比較

• 因子分析は観測された変数間の共分散と相関を少数の共通因子で説明．主成分分析は

観測された変数の分散を説明

• 因子分析で因子の数を変えると他の因子も変わる．

主成分分析で主成分の数が変わっても，最初のほうの主成分は変わらない

• 因子分析は複雑で多義的．主成分分析は明確な結論が出る．

(76)

数量化

I

類

次のようなデータがある．

X1

X2

Y

Med

G2

9.58 Med

G3

6.29 Med

G4

12.35 Hi

G1

9.12 Hi

G2

13.84 Lo

G3

1.63 Hi

G4

15.37 Hi

G2

13.45

8 つの「物件」があり，それぞれが

2 つのランキング

X1, X2

によって分類されている．

X1

の方では

Lo, Med, Hi

の

3 段階，

X2

の方では

G1

∼

G4

の

4 段階になっている．そ

して，

Y

はこれらの物件の価格というふうに考えることにしよう．

(77)

たとえば，一番上の物件は，

Med, G2

というランキングで，価格は

9.58 である．

数量化

I

類とは

この場合，

X1, X2

のランキングは数値ではなく，

ファクター

である．因子とかカテゴ

リカル変数

,

あるいは名義変数ともいう．

このように，ファクターを説明変数としてもつデータがあった場合，それらを数量化し

て，重回帰分析で扱おうというのが数量化

I

類である．

R

には，数量化

I

類そのものの関数は実装されていない．しかし，すでに出てきた重回

帰分析をファクターについて適用することで，数学的には同等の分析ができる．

(78)

■傾向を目で読み取ってみる

データをよく見ると，一定の傾向が見えてくる．ここで

X = Med

に固定してみると，

G2 < G3 < G4

という序列であるし，逆に見ていくと，

Lo < Med < Hi

という序列が

ある．

そこでたとえば，

Lo G4

だったら，

Y

はどうなるかを予測することが可能なはずであ

る．それが数量化

I

類の役割だ．

(79)

R

による処理

以下のように，ごく短いスクリプトで処理が完了する．

D <- read.table("NumerizeI.txt",header=T)

#

データ読み込み

## Y

を目的変数，

X1

と

X2

を説明変数として線形モデル

(lm)

を適用

result <- lm(Y ~ X1 + X2, data = D)

summary(result)

#

結果のまとめを出力

(80)

結果の出力

結果は重回帰分析のものと同じである．ただし，

predict

関数で予測値を出すと，一般

的な数量化

I

類と同じ結果になっている．

■summary の出力

Residuals:

1

2

3

4

5 -2.986e-01

2.776e-17

2.986e-01

1.665e-16

3.443e-01

6

7

8 2.776e-17 -2.986e-01

-4.571e-02

1

2

3

4

5 -2.986e-01

2.776e-17

2.986e-01

1.665e-16

3.443e-01

6

7

8

(81)

Coefficients:

Estimate Std. Error t value Pr(>|t|)

(Intercept)

9.1200

0.4405

20.704 0.00232 **

X1Lo

-8.2771

0.7446 -11.117

0.00799 **

X1Med

-3.6171

0.4078

-8.870

0.01247 *

X2G2

4.3757

0.5265

8.311 0.01417 *

X2G3

0.7871

0.7446

1.057 0.40126

X2G4

6.5486

0.5767

11.355 0.00767 **

---Signif. codes:

0 ‘

***

’

0.001 ‘

**

’

0.01 ‘

*

’

0.05 ‘

.

’

0.1 ‘ ’

1 Residual standard error: 0.4405 on 2 degrees of freedom

Multiple R-squared:

0.9973,Adjusted R-squared:

0.9907

(82)

■predict の出力

1

2

3

4

5 9.878571

6.290000 12.051429

9.120000 13.495714

6

7

8 1.630000 15.668571 13.495714

この出力を与えたデータと照らし合わせると，

Med G2

となっている物件では，価格が

9.88 という予測値になるということがわかる（以下同様）．

(83)

クラスター分析

多数の「もの」に対して数量的なデータを元にした分類を行うことは有

益な情報処理の手段である．

ものとものの間の「距離」を使うことによって，近いもの＝類縁という

関係を抽出して，グルーピングを行う手法がクラスター分析である．結

果を階層的なツリー構造で表せる階層的クラスター分析と，非階層的ク

ラスター分析とがある．

ここでは視覚的にわかりやすい階層的クラス

ター分析を紹介する．

(84)

距離を定義する

平面上の

2 点

(x

1

, y

1

), (x

2

, y

2

)

があったとすると，ピタゴラスの定理から

d =

√

(x

1

− x

2

)

2

+ (y

1

− y

2

)

2

三次元空間の

2 点

(x

1

, y

1

, z

1

), (x

2

, y

2

, z

2

)

ならば，

d =

√

(x

1

− x

2

)

2

+ (y

1

− y

2

)

2

+ (z

1

− z

2

)

2

n

次元に拡張し，座標

(

変量

)

のベクトル

x

1

, x

2

, . . . x

p

があるとした場合，

x

1

, x

2

の距

離は

d =

√

(x

1,1

− x

1,2

)

2

+ (x

2,1

− x

2,2

)

2

+

· · · + (x

n,1

− x

n,2

)

2

のようになる．多次元の長さというのは想像するのが難しいが，なんとなく低次元のイ

メージで押さえておけばよい．

このような定義による距離を

ユークリッド距離

という．

(85)

なお，こんなふうに書くとややこしい．既に出てきたノルムを使って，

d =

||x

1

− x

1

||

このように表現する．

要するに，

2 つのデータの間の「近さ」を，幾何学的な発想で定義している

(

他の定義も

ある

)

．

クラスター分析では，データの間の「近さ」を定義して，

それによって近親関係をまとめていく．

(86)

クラスター分析の実例：各年齢別の平均余命の国別比較

m0

m25

m50

m75

w0

w25

w50

w75

Algeria

63

51

30

13

67

54

34

15 Cameroon

34

29

13

5

38

32

17

6 Madagascar

38

30

17

7

38

34

20

7 Mauritius

59

42

20

6

64

46

25

8 Reunion

56

38

18

7

62

46

25

10 Seychelles

62

44

24

7

69

50

28

14 South Africa(B)

50

39

20

7

55

43

23

8 South Africa(W)

65

44

22

7

72

50

27

9 Tunisia

56

46

24

11

63

54

33

19 Canada

69

47

24

8

75

53

29

10 Costa Rica

65

48

26

9

68

50

27

10 Dominican Rep

64

50

28

11

66

51

29

11 . . .

. . .

Ecuador

57

46

28

9

60

49

28

11 全

31 ヵ国．

mxx, wxx

は

xx

歳の男性

(

女性

)

の平均余命

(87)

R

による処理

##

データフレームに読みこむ

life <- read.table("LifeExp.txt",header=T)

##

行の名前（表の左端）を

country

に代入

country <- row.names(life)

##

各国のデータの間のユークリッド距離を計算

dist <- dist(life)

postscript("LifeExp.ps",horizontal=F, width=7,

height=7,onefile=TRUE)

##

最長距離法でクラスター分析

plot(hclust(dist, method = "complete"),

labels = country,

#

データごとのラベルを国名にする

(88)

得られた結果

Trinidad(62) Canada United States (W66) Argentina South Africa(W)

United States (66) United States (67)

Trinidad (67)

United States (NW66)

Chile

Seychelles Grenada Jamaica

Reunion

Mexico

Colombia Honduras Mauritius _Greenland

Algeria

Costa Rica Panama

Dominican Rep Nicaragua Tunisia El Salvador Ecuador South Africa(B)

Guatemala Cameroon Madagascar

0 10 20 30 40 50 60 最長距離法で得られたクラスターツリー hclust (*, "complete")国名距離

(89)

クラスター化のアルゴリズムの比較

Cameroon Madagascar Trinidad(62) Tunisia South Africa(B) Guatemala Algeria Reunion Seychelles Costa Rica Panama Argentina South Africa(W)

Canada United States (W66) Dominican Rep Nicaragua United States (NW66) Honduras Chile Trinidad (67)Grenada Jamaica

Ecuador

El Salvador Mexico Colombia Mauritius Greenland

0 5 10 15 20 (a) 最短距離法 hclust (*, "single")国名距離 Trinidad(62) Canada United States (W66) Argentina South Africa(W)

Trinidad (67) United States (NW66) Chile Seychelles Grenada Jamaica Reunion Mexico Colombia Honduras Mauritius Greenland Algeria Costa Rica Panama Dominican Rep Nicaragua Tunisia El Salvador Ecuador South Africa(B)

Guatemala Cameroon Madagascar

0 10 20 30 40 50 60 (b) 最長距離法 hclust (*, "complete")国名距離 Cameroon Madagascar South Africa(B) Guatemala Trinidad(62) Algeria Tunisia

Canada United States (W66) South Africa(W) Argentina Costa Rica Panama Dominican Rep Nicaragua Reunion El Salvador Ecuador Honduras Mexico

ColombiaMauritius Greenland

Seychelles

Chile

Grenada Jamaica Trinidad (67)

United States (NW66) 0 10 20 30 40 (c) 群平均法 hclust (*, "average")国名距離

(90)

R R-console R R Rscript R-console GUI 1

多変量解析入門

R

超入門

インタラクティブなインターフェース

R-console

の画面を呼び出してコマンドを打ち込む

ターミナルから

R

と打ち込むと，上と同様のキャラクター画面が現れる

ので，コマンドを打ち込む

バッチ処理のためのインターフェース

エディタで

R

のソースファイルを作って，

Rscript

で読み込んで実行

一見むつかしそうだが，

R-console

や ターミナルで打ち込むコマンドを

並べるだけなので，この方式に慣れる方がよい．保存性や編集の容易さ

などで計り知れないメリットがある．プログラミングの経験のある人な

ら，この方式を最初から取るほうがよい．

また，

GUI

のメニューで操作するよりも．

1

行ずつ処理を書いていくとい

うやり方のほうが，手続きを言語化することができる分，分かりやすい．

これは大きな利点．

2

組のデータの相関から

２変量のデータ

説明変数：

x = [x

, x

, . . . , x

]

目的変数：

y = [y

, y

, . . . , y

]

太字はベクトルであることを表す．つまり，これらのデータは

n

次元の

ベクトル．多変量解析においては．説明変数が

x

, x

, . . .

となる．

行列との積を作る場合には，しばしば縦ベクトルとして扱われる．

x =

x

x

...

x

, y =

y

y

...

y

スペースがもったいないので転置を意味する

T

を付けて次のように書く

こともある．

x = [x

, x

, . . . , x

]

2

変量から作られる統計量

x, y

の２系列のデータから，２つの平均，２つの分散（標準偏差），そし

平均

¯

x =

1

n

(x

やターミナルで打ち込むコマンドを