六本木で働くデータサイエンティストのブログ

元祖「銀座で働くデータサイエンティスト」です / 道玄坂→銀座→東京→六本木

なぜ項目ごとに単純な集計をするより、多変量解析(重回帰分析)をした方が正確な結果を返すのか

ちょっと前の記事(単純な集計とデータサイエンスによる分析とで結果が食い違うかもしれない3ケース)に裏先生からツッコミを頂き、慌てて学部1年生の頃の教科書を開いて復習しまして。。。いやー、忘れてることが多過ぎて大変(汗)。知るは一時の恥というのをつくづく思い知りましたとさ。本当に裏先生ありがとうございました。


ということで、その復習内容の確認と同時に、あの時の裏先生のツッコミ内容をかみくだいて紹介するのも兼ねて、ここはひとつざっくり書いてみようかと思います。


項目ごとの単純集計は「単変量」解析(もっと言えば単相関)に過ぎず、多変量データ全体のことは分からない


前回用いたサンプルデータは、基本的にはa1-a7が0or1のみの二値で構成される事実上のカテゴリカルデータ*1で、cvも"Yes"or"No"のみの二値で構成されるカテゴリカルデータです。


で、二値のカテゴリカルデータだけで構成されているので、これらを例えば"Yes"/"No"ごとにa1-a7でまとめて平均値を出して差を取るとか、はたまた"a6"/"a7"ごとに"Yes"/"No"でまとめて平均値を出して差を取るという操作は、そのまま相関係数を求めるのとほぼ同じ結果になります。前者の「平均値を出して」「差を取る」というのはSQLでクエリを書くだけで簡単に実装できる計算なので、どこのデータ分析の現場でも広く用いられているやり方ではないかと思います。


実際にcvの列を0or1のdouble型に直してcv2とし、例えばxmというマトリクスに改めてその上で確認してみると、

a1 a2 a3 a4 a5 a6 a7 CV
0.401 0.583 0.479 0.942 0.307 0.056 0.500 No
0.605 0.417 0.494 0.436 0.684 0.927 0.493 Yes
0.203 -0.166 0.015 -0.506 0.377 0.871 -0.007 Yes - No
> cor(xm)
              a1           a2           a3           a4           a5          a6           a7          cv2
a1   1.000000000 -0.060005093 -0.003174519 -0.114786423  0.100728867  0.16279221  0.006040243  0.203336993
a2  -0.060005093  1.000000000  0.002018496  0.100540928 -0.069996842 -0.15067652 -0.027329502 -0.166000037
a3  -0.003174519  0.002018496  1.000000000  0.007289939 -0.017589243  0.02156622 -0.030856103  0.014671884
a4  -0.114786423  0.100540928  0.007289939  1.000000000 -0.195663709 -0.48436622  0.009203000 -0.546551028
a5   0.100728867 -0.069996842 -0.017589243 -0.195663709  1.000000000  0.33523894 -0.004062489  0.377349769
a6   0.162792208 -0.150676518  0.021566219 -0.484366216  0.335238942  1.00000000 -0.011446289  0.871454377
a7   0.006040243 -0.027329502 -0.030856103  0.009203000 -0.004062489 -0.01144629  1.000000000 -0.006666815
cv2  0.203336993 -0.166000037  0.014671884 -0.546551028  0.377349769  0.87145438 -0.006666815  1.000000000 # この行を見て下さい


cv2の行と見比べてみると、a4だけちょっとずれてますが基本的には合っているはずです。同じことは、例えばExcelで「分析ツール→相関or共分散」もしくはCORREL関数でできます。


しかし、これはあくまでも「単相関」でしかなく、cv2という単一の水準から見た値でしかありません。ここに注意が必要です。


もっと書いてしまうと、今回の例で言えばこのやり方だとcv2から見た「a1-a7からの見かけの影響」しか見えないことになります。すなわち、a6とa7からcv2への影響が実は互いにだだ被りしていたとしても、これを見分ける術はないわけです。


これが、SQLで簡単に実装できるような項目ごとの単純集計をしているだけでは、正確な結果にならない可能性があるということの理由です(この2つ後の「偏相関について」の項もお読み下さい)。


(多変量解析の一例として)重回帰分析は項目間の影響を排除した「偏」回帰係数を出すので、より正確に項目ごとの影響度を量れる


せっかくなので、学部1年生の頃の教科書に立ち返ってみましょう。僕がその頃使っていたのは、古典的名著『統計学入門』です。この記事書く前に見たら、僕の字でシャーペンか何かで書き込んだ跡がありましたw


統計学入門 (基礎統計学)

統計学入門 (基礎統計学)


さて本題ですが。大事なポイントとして、そもそもの大前提&問題意識として「多変量データはそのまま扱うと説明変数(つまり項目)同士が影響し合ってしまう」「そのまま単変量のつもりで分析すると結果が歪む」という上に書いた通りの話があり、「それを回避するためのモデルを予め仮定しこれを解く」というその後の流れが想定されています。


言い換えると、「多変量になった時点で単純集計するとヤバいかも」という危機意識があるということですね。


なので、モデルの立て方も「各説明変数の被説明変数への(他の説明変数の影響を除いた純粋の)影響度を線形結合したモデル」という形を取っています。予めモデル自体に「説明変数同士に生じ得る相互影響を排除しておく」という意味合いがある、ということです。こうしておけば、そのモデルを単純に数学的に正確に求める*2ことでそもそも説明変数同士の相互影響が一切ないモデルが得られ、結果的に正しい結論に至るわけです。


その方法論は『統計学入門』pp.271-272にも載っていますが、web上でも閲覧できるものとしては有名な青木繁伸先生(群馬大)のサイト『おしゃべりな部屋』の「統計学自習ノート」の偏回帰係数の求め方に関するページを引用しておきます。ここでのやり方は、原理的には通常の線形重回帰分析でもロジスティック重回帰分析でも変わりません*3


要は、いわゆる最小二乗法(Ordinary Least Square: OLS)に則って、1) 被説明変数を説明変数と偏回帰係数の線形結合から成る回帰式を立て、そこから2) 誤差項の式に変え、3) 誤差項の平方和をそれぞれのn個の説明変数に対して偏微分し、4) 出来上がったn本の連立方程式をn個の偏回帰係数について解く、というやり方です。


これを『統計学入門』pp.271-272に従って書くと、

1) 回帰式を立てる

Y_i = \beta_1 + \beta_2 X_{2i} + \beta_3 X_{3i} + \cdots + \beta_k X_{ki} + \eps_i

2) 誤差項の式に変える

\eps_i = Y_i - (\beta_1 + \beta_2 X_{2i} + \beta_3 X_{3i} + \cdots + \beta_k X_{ki})

3) 誤差項の平方和を個々の説明変数で偏微分して0と置く

2)の誤差項の平方和

S = \Sigma \eps_i^2

に対して、その説明変数による偏微分を0と置く。すなわち

\frac{\partial S}{\partial \beta_1}=0, \frac{\partial S}{\partial \beta_2}=0, \cdots, \frac{\partial S}{\partial \beta_k}=0

4) 偏微分して得られた連立方程式を解く

これはごく普通のn元1次連立方程式なので、単純に解ける。そこで得られた

 \hat{\beta}_1, \hat{\beta}_2, \cdots , \hat{\beta}_k

こそが求める偏回帰係数となる。


この計算自体は、式変形の結果として行列計算に直すことができます。実際に行列計算として数値的に解いている分析ツールは多いです。


またOLSより進んだ方法として、最尤法(Maximum Likelihood Estimation)で偏回帰係数を求めることもできます。ここでは簡単のために割愛しますが、最尤法で重回帰分析を行う方法についてはid:uncorrelatedさんの記事(Rと手作業で覚える最尤法)が分かりやすくて良いです。glm(){base}関数は最尤法を用いて偏回帰係数を計算しています。


こうして得られた偏回帰係数 \hat{\beta}_1, \hat{\beta}_2, \cdots , \hat{\beta}_kは、実際に上の方で事前に定めた「各説明変数の被説明変数への(他の説明変数の影響を除いた純粋の)影響度を線形結合したモデル」という前提を満たして求められています*4


なので、重回帰分析(もっと一般的に言うなら多変量解析)で求めた結果の方が、単純集計とは異なり説明変数同士の影響を混同することなく、正確な値を返すというわけです。


Rでの計算例は前回の記事でも紹介した通りですが、一応再掲しておきます。

# サンプルデータは"sample_d"というデータフレームに入れてあるものとする

> sample_d.glm<-glm(cv~.,sample_d,family="binomial")
> summary(sample_d.glm)

Call:
glm(formula = cv ~ ., family = "binomial", data = sample_d)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-3.6404  -0.2242  -0.0358   0.2162   3.1418  

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept) -1.37793    0.25979  -5.304 1.13e-07 ***
a1           1.05846    0.17344   6.103 1.04e-09 ***
a2          -0.54914    0.16752  -3.278  0.00105 ** 
a3           0.12035    0.16803   0.716  0.47386    
a4          -3.00110    0.21653 -13.860  < 2e-16 ***
a5           1.53098    0.17349   8.824  < 2e-16 ***
a6           5.33547    0.19191  27.802  < 2e-16 ***
a7           0.07811    0.16725   0.467  0.64048 # ←コレだよコレ!!!    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 4158.9  on 2999  degrees of freedom
Residual deviance: 1044.4  on 2992  degrees of freedom
AIC: 1060.4

Number of Fisher Scoring iterations: 7


なお、簡潔な説明が『おしゃべりな部屋』の重回帰分析の「多重共線性」パートのところにあります*5。人によってはこちらの方が分かりやすいかも。


当然の話ですが、いかに多変量解析を使おうとも無視できないほど説明変数の多重共線性が強ければ、偏回帰係数の値は歪んでしまいます。説明変数に可能な限り「互いに独立」なものだけに限るというのは、分析以前のポイントとして肝に銘じておく必要があります。


ちなみに、線形モデル&最小二乗法周りを再確認したい人は、『統計学入門』の続編に当たる『自然科学の統計学』あたりがお薦めです。


自然科学の統計学 (基礎統計学)

自然科学の統計学 (基礎統計学)


上記の偏回帰係数を求めるくだりが、『自然科学の統計学』pp.36-40に出てきます(連立偏微分方程式から行列計算に落としていく導出過程も全部書いてあります)。


今回取り上げた本はどちらもかなり古い部類に属しますが、Rを用いたデータ分析の現場で使われる関数のかなりの部分の統計学的基礎を網羅しているので、脇に置いておいて決して損はしないと思います。


偏相関について


裏先生からのツッコミにもありましたし、『統計学入門』pp.52-54にも書いてありますが、『おしゃべりな部屋』の「偏相関係数」のページの説明が手っ取り早いでしょう。


要は、一番最初に示した「単純集計の結果がおかしい」理由は、a6とa7とが無視できない相関を持っていて、そのせいでcv2へのa6とa7からの影響が混ざってしまったからです。


それをシンプルに解決して、相関係数を求めた時と同じように計算する方法が偏相関です。Rなら{ppcor}パッケージのpcor()関数で計算できます。

> pcor(xm)
$estimate
              a1           a2           a3           a4           a5          a6           a7          cv2
a1   1.000000000 -0.027097667 -0.004665751 -0.005058604  0.026653584 -0.03065927  0.006356354  0.111647484
a2  -0.027097667  1.000000000  0.003302999  0.011735712 -0.007277052 -0.01330482 -0.028772727 -0.058979223
a3  -0.004665751  0.003302999  1.000000000  0.019140545 -0.025398517  0.01818792 -0.030596755  0.002985735
a4  -0.005058604  0.011735712  0.019140545  1.000000000  0.014636602 -0.02011215  0.007394668 -0.284084473
a5   0.026653584 -0.007277052 -0.025398517  0.014636602  1.000000000  0.01549824 -0.002781801  0.175562475
a6  -0.030659270 -0.013304820  0.018187916 -0.020112145  0.015498244  1.00000000 -0.010926019  0.807244783
a7   0.006356354 -0.028772727 -0.030596755  0.007394668 -0.002781801 -0.01092602  1.000000000  0.005953911
cv2  0.111647484 -0.058979223  0.002985735 -0.284084473  0.175562475  0.80724478  0.005953911  1.000000000


Excelにはデフォルトでは入っていないので、上記『おしゃべりな部屋』のページを参照するなどして自前で実装する必要がありそうです。


細かいことですが


裏先生からのご指摘2点ほどについてお答えをば。ちなみにどちらも統計学的なポイントというより、現実のビジネス的な(特にweb業界での)データサイエンスにおけるポイントと言うべきかもですが。

単純な分析で総合的な分析(多変量解析)の結果は予測できないということ。


「何故かある」のではなく,「必然的にある」。ごくたまに単純な分析結果から予想した結果が総合的な結果と一致することもあるのだが。


概して、web業界でこの手の影響要因の分析を行う際に説明変数として選ばれるものは、サービス運用上の経験や前提知識などに基づいて、互いの独立性が高いケースが多いです。このため、単純集計で求めた結果でもそれなりに正確だったりします。


実際、「単純な分析結果から予想した結果が総合的な結果と一致すること」は、統計学的にはたまたまでしかないとは言っても、データ分析の現場ではむしろ多いという印象です。時々僕らのような人間が念のために多変量解析をかけたら、食い違いが判明するというケースがまばらにある程度です。


とは言え、これを読んでいるデータ分析実務者の諸氏で、普段は単純集計しか使わないという方々には、「単純集計しかやっていないといつかどこかで誤った結果に惑わされかねない」ということを是非肝に銘じてもらえればと思います。単純集計だけならSQLで簡単に自動化して実装できるとはいえ、時々はローカルに生データを落としてきて偏相関ぐらい見てもバチは当たらないということで。


なお、似たような話になりますが、

「a7はa6との組み合わせにおいてのみCVR増に貢献していると考えられるので、a7はa6との間にどういう関係性があるのかをもっと細かく調べるべき」と答えることになると思います」と述べている。


わかっているのかどうかよくわからない。


これは「今回のような結果が出た時にデータ分析担当者がビジネスサイドでの意思決定者に返すべきコメント」の一例であるに過ぎません。もちろん、データ分析担当者の側ではこれがa6とa7との間の(一部に限って)相関が高いせいであると分かっているわけですが、それ自体が分かってもビジネス上は何の役にも立たないわけです。


大事なことは、こういう問題を見つけたら「a6とa7との間には何か見落としている重大な関連性があるのではないか」「a6とa7は同時に分析指標には入れられない」「a6とa7のどちらを残すのか」「他にもっと独立性の高い役立つ指標はないのか」とデータ分析担当者がビジネスサイドに問うこと、でしょう。それは統計学上の議論の向こう側にある話だと理解しています。


最後に


基礎の基礎ほど習った時期が古過ぎる上に、最近になればなるほどあれやこれやの手法をツールやライブラリを使ったりして何気なく使ってしまうので、忘れやすいんですよね。。。時々昔の教科書開いて、復習するのってものすごーーーく重要だなとつくづく思いました。これからも基礎を忘れないよう、暇を見て復習せねば。

*1:厳密には「名義尺度」と呼ぶべきなんでしょうが

*2:計算プロセスで説明変数同士が混同しないようにするということ

*3:ロジスティック回帰は基本的には通常の回帰分析を変数変換したものに過ぎないので

*4:偏微分によってある程度独立に個々の説明変数ごとの影響度を表現しているという側面もある気がする

*5:ちなみに前回のケースで言えば、a6-a7の間の偏相関は実際には共線性を疑わせるほど強くなく設定してあったり。。。