静粛に、只今統計勉強中

仕事でデータ分析をすることになったバリバリ文系アラフィフのおっさんが、独学で統計の勉強を始めました。

回帰分析における対数変換の意味を実感してみよう1

2018年11月試験で統計検定2級に合格したので、今度は準1級に挑戦しようと思い、 

日本統計学会公式認定 統計検定 2級 公式問題集[2012~2014年]

日本統計学会公式認定 統計検定 2級 公式問題集[2012~2014年]

に掲載されている<例題>を解いてみました。

準1級は、2015年6月から始まったので、上の問題集が出版された時点では未開始でした。なので、<例題>は「こんな感じの問題を出しますから、参考にしてくださいね」という日本統計学会の配慮といったところでしょうか。

まだ、実際に出題された過去問は見ていませんが、<例題>は、ほとんどが2級レベルでした。こんなことなら、2級を受ける前にも解いておけばよかった。

<例題>で2級のレベルを超える問題は36問中10問(うち論述が2問)ありましたが、全然全く歯が立ちませんでした。残り10カ月で受験できるようになるか、かなり不安です。

 

それはさておき、今回のお題は回帰分析における対数変換です。というのも、<例題>の問11で取り上げられたデータが、その意味を実感するのにうってつけと思えたから。

取り上げられたのは、映画館に関する経済産業省「平成22年特定サービス産業実態調査(確報)」のデータです。*1

散布図にすると、こんな感じ。

f:id:cyclo-commuter:20190828144616p:plain

横軸が映画館の従業者数、縦軸がスクリーン数です。

決定係数が0.9198と極めて高いことから、この回帰モデルは有効のように見えます。

しかし、問11〔2〕の問題は、散布図と回帰診断図から判断して問題点を述べよ、というもの。回帰分析において、決定係数が高くてもダメな場合があるなんて、2級試験の問題では出てきませんでしたよねえ。

 

ともあれ、データもあることですし、Rで分析の流れを追体験してみましょう。

 

まずは、経済産業省の『「映画館」統計表一覧データ』の都道府県別第10表から、従業者数とスクリーン数のデータを抜き出します。

さらに、余計な行を削除します。
島根県徳島県のスクリーン数が秘匿されていますので、これも削除。
あと、カンマがあるとRで文字列扱いされちゃうので、数値の書式を標準にしておきます。

f:id:cyclo-commuter:20190828150101p:plain

Range(A1:C46)をコピーしておきます。

 

Rコマンダーを起動して、[データ]→[データのインポート]→[テキストファイルまたはクリップボード,URLから…]をクリック。

f:id:cyclo-commuter:20190828145340p:plain

クリップボード」にチェックを入れて、【OK】をクリック。

f:id:cyclo-commuter:20190828145702p:plain

[統計量]→[モデルへの適合]→[線形回帰]をクリック。

f:id:cyclo-commuter:20190828150607p:plain

目的変数にスクリーン数、説明変数に従業者数を選択し、【OK】をクリック。

f:id:cyclo-commuter:20190828150733p:plain

計算結果がこちら。

Call:
lm(formula = スクリーン数 ~ 従業者数, data = Dataset)

Residuals:
    Min      1Q  Median      3Q     Max 
-54.566  -9.848  -4.421   6.684  82.728 

Coefficients:
             Estimate Std. Error t value  Pr(>|t|)    
(Intercept) 18.465565   3.862580   4.781 0.0000207 ***
従業者数     0.106381   0.004791  22.204   < 2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 20.41 on 43 degrees of freedom
Multiple R-squared:  0.9198,	Adjusted R-squared:  0.9179 
F-statistic:   493 on 1 and 43 DF,  p-value: < 2.2e-16

<例題>の問11と同じ値になりました。

続いて、回帰診断図を出力します。

 

[モデル]→[グラフ]→[基本的診断プロット]をクリック。

f:id:cyclo-commuter:20190828151532p:plain

結果がこちら。

f:id:cyclo-commuter:20190828151740j:plain

 これらの図で何を診断するかというと、回帰モデルの誤差\varepsilonが持つべきとされる以下の性質です。

  1. 不偏性:それぞれの誤差の期待値は0である( E(\varepsilon)=0
  2. 等分散性:それぞれの誤差の分散は一定値である(V(\varepsilon)=\sigma^2
  3. 独立性:互いに異なる誤差は独立である( Cov(\varepsilon_i,\varepsilon_j)=0,i\neq j
  4. 正規性:誤差は正規分布に従う

『RとRコマンダーではじめる多変量解析』荒木孝治編著(日科技連出版社)P22

 これらの性質を満たさない場合、元のデータにバイアス(偏り)があったりして、得られた分析結果が信頼できないものになってしまうのだそうです。

上図の場合、

  • Residuals vs Fitted は、もっと満遍なくプロットされるはずが、予測値が大きくなるほど残差の絶対値が大きくなっています。(等分散性に問題あり)
  • Normal Q-Q は、直線状にプロットされるはずが、東京・福岡・愛知が直線から大きくはみ出しています。(正規性に問題あり)
  • Scale-Location は、もっと水平にプロットされるはずが、予測値が大きくなるほど標準化残差の絶対値の平方根が大きくなっています。(等分散性に問題あり)
  • Residuals vs Leverage は、クックの距離0.5の線の内側に収まるべきが、愛知が0.5以上、東京が1以上になっています。(東京と愛知は外れ値っぽい)

ということなんですね。うーむ。

散布図ではそこそこ回帰直線に近いところに集まっているし、決定係数も0.9を超えるくらい高い。けれども、回帰モデルの信頼性は低い。ということがあるんですねえ。

 

<例題>の問11の解答はここまでなんですが、じゃあどうすれば信頼できる回帰分析になるんでしょうか?