Last updated: October 2009. Kajiyama                    [ 目次に戻る ]

三種類の触媒で作られた製品の製造原価比較: 一元配置の分散分析と多重比較

(対応なし・標本数が異なる)


  1. 「対応なし・標本数が異なる」場合の課題

    三以上の群の平均値を比較するのであるが,各群の標本数が同じでない場合がある.

    例えば実験に失敗してデータを取れない場合や,予定した被験者が当日いなかった場合など.実験装置をゴジラが壊したなど.

    このような,「各群の標本数が異なる場合の分散分析」をあなたは大切にする.



    三種類の触媒で作られた製品の製造原価を繰り返し調べたものが表1である.

    「三つの種類の触媒(水準A1,水準A2,水準A3)の製造原価(円)の母集団の平均値に差があるといえるか」,有意水準5%で検討する.


    表1 三種類の触媒で作られた製品の製造原価 「標本数が異なる アンスタック形式」

    工業関係の分散分析法の解説は,「分散分析法入門,石川馨,米山高範,日科技連出版社,1967.」が詳しい.表1のデータはp39.


    サンプル
    01
    02
    03
    04
    05
    06
    07
    08
    平均値
    触媒・水準A1 触媒・水準A2 触媒・水準A3
    16214 15997 17018
    16785 16254 16325
    16009 15432 16485
    17532 15886   
    16532      
    17489      
    15963      
    16986      
    16680 15890 16610


  2. 分散分析の関数は「スタック形式」のデータで処理する.


  3. Rに「スタック形式」のデータを読み込み,データの確認を行う


    表2 三種類の触媒で作られた製品の製造原価 「標本数が異なる スタック形式」

    サンプル
    01
    02
    03
    04
    05
    06
    07
    08
    09
    10
    11
    12
    13
    14
    15
    製造原価 触媒
    16214 水準A1
    16785 水準A1
    16009 水準A1
    17532 水準A1
    16532 水準A1
    17489 水準A1
    15863 水準A1
    16986 水準A1
    15997 水準A2
    16254 水準A2
    15432 水準A2
    15886 水準A2
    17018 水準A3
    16325 水準A3
    16485 水準A3


    (1) 「スタック形式」での操作の手順

    > y <- read.table ("clipboard", header=TRUE )
    > y
       製造原価   触媒
    1     16214  水準A1
    2     16785  水準A1
    3     16009  水準A1
    4     17532  水準A1
    5     16532  水準A1
    6     17489  水準A1
    7     15863  水準A1
    8     16986  水準A1
    9     15997  水準A2
    10    16254  水準A2
    11    15432  水準A2
    12    15886  水準A2
    13    17018  水準A3
    14    16325  水準A3
    15    16485  水準A3
    > str (y)
    'data.frame':   15 obs. of  2 variables:
     $ 製造原価: int  16214 16785 16009 17532 16532 17489 15863 16986 15997 16254 ...
     $ 触媒    : Factor w/ 3 levels "水準A1","水準A2",..: 1 1 1 1 1 1 1 1 2 2 ...
    


    (2) データフレームの読み込み

    ここでは,表シートあるいはエクセルの表計算シートからクリップボード経由で読み込む操作を行う.

    表2 三つの種類の触媒で作られた製品の製造原価の緑の部分をマウスで選択し,これをコピーする.

    次に,Rの「コンソール」画面に,

    y <- read.table ("clipboard", header=TRUE

    と記述し,コピーした[clipboard]データファイルを,Rの内部のデータフレーム,ファイル名x に直接読み込む.

    header=TRUE は第1行が列の変数名になっていることを指示している.

    ) #    データをコピーした後で,カッコを入力する.ENTERーキーを押す.


    (3) 読み込んだデータフレームの確認

    y     と記述し,Rの内部に作成したデータフレーム y を,表示し確認する.


    (4) データフレームの構造確認

    str (y)  と記述し,Rの内部のデータフレーム y の2変数の内容を表示確認する.

    Rが理解した,製造原価と触媒の2変数名のリストと変数ごとのカテゴリの数とカテゴリ名が示される.

    ・分析で使用する変数名とカテゴリ名はRが理解した名前を使用する.


    注意: 水準A1〜A4を含む因子(触媒)のデータ形式が [Factor] 形式でないときは,正しい統計処理ができなくなる.

    y$触媒 <- factor( y$触媒 )  と記述し,データフレームの水準の変数を[Facter]形式に変更する.



  4. 1因子Aの基礎統計


    a. 因子Aの基礎統計量

    Rの「コンソール」画面に,

    by ( y$製造原価, y$触媒, summary)   と記述し,「スタック形式」のデータで,群別の基礎統計利用を計算する.

    > by ( y$製造原価, y$触媒, summary)
    y$触媒: 水準A1
       Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
      15860   16160   16660   16680   17110   17530 
    ------------------------------------------------------------------------ 
    y$触媒: 水準A2
       Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
      15430   15770   15940   15890   16060   16250 
    ------------------------------------------------------------------------ 
    y$触媒: 水準A3
       Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
      16320   16400   16480   16610   16750   17020 
    


    b. 因子Aの平均値の推移図    # plotMeans関数を使用するにはパッケージ 「Rcmdr」を読み込んでおくこと

    boxplot ( y$製造原価 ~ y$触媒)  と記述し,箱ひげグラフを表示する.

    plotMeans ( y$製造原価, y$触媒, error.bars="se")  と記述し,平均値の推移図を表示する.
     
    c. グラフの観察



  5. バートレット(Bartlett)の検定--分散の等質性(等分散性)の検定


    3水準(3群)の分散が等しいか,あるいは等しくないかを検定する.3水準のデータの分散は等質であるとの帰無仮説を立てる.

    分散分析は,分散が等しいという仮定の上で構築されているので標本の測定値にもとづいて分散が一様であるか検討する.

    ただし,分散分析のF検定はその前提条件が崩れても,検定結果は信頼できるといわれている.


    Rの「コンソール」画面に,

    bartlett.test ( y$製造原価 ~ y$触媒 ) # と記述し,「スタック形式」のデータで,分散の等質性の検定を行う.

    バートレットの検定は水準ごとの標本の数が同じでない場合も使用できる.

    p値が0.05以下のときに水準の分散は等しくないと結論する.


    >  bartlett.test ( y$製造原価 ~ y$触媒 )
    
            Bartlett test of homogeneity of variances
    
    data:  y$製造原価 by y$触媒 
    Bartlett's K-squared = 1.6115, df = 2, p-value = 0.4468
    

    分析結果

    p値が0.4468であり,p値>0.05であるので各水準の分散は等しいと結論し,分散分析に進む.



  6. 「対応なし・標本数が異なる」場合の一元配置・分散分析


    「三種類の触媒(A1〜A4)で製造した製品の製造原価の母集団の平均値に差があるといえるか」,有意水準5%で検討する.

    帰無仮説は「各群の母集団の製造原価の平均値は等しい」.

    対立仮説は,「全ての群の組み合わせに平均値の差がある」のでなく,「各群の少なくとも一組に平均値の差がある」.



    Rの「コンソール」画面に,

    summary ( aov ( y$製造原価 ~ y$触媒 )) # と記述し,「スタック形式」のデータで,分散分析を行う.

    summary ( aov ( 製造原価 ~ 触媒, y )) # の書式も同じ結果となる.

    > summary ( aov ( y$製造原価 ~ y$触媒 ))
    
                Df    Sum Sq  Mean Sq  F value  Pr(>F)  
    y$触媒       2   1729797  864899   3.0007   0.08775 .
    Residuals   12   3458741  288228                  
    ---
    Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 
    有意性の記号 :    ***  (0.1%で有意)  **  (1%で有意)   *  (5%で有意)   . (有意の傾向がある)
    
    
    分散分析表が表示される.これを整形すると,以下のようになる.
    -----------------------------------------------------------          (変動)  (不偏分散) (分散比) 変動要因 自由度 平方和  平均平方和 F値 P値 ----------------------------------------------------------- 因子A 2    1729797 864899  3.0007  0.08775 誤差E 12    3458741 288228 -----------------------------------------------------------
    分析結果

    因子(触媒)AのA1〜A3の水準間の変動について,P値>0.05であることから,5%水準で帰無仮説は棄却されない.

    因子Aは有意差がない.製造原価の平均値は触媒の種類間に違いがあるとはいえない.



  7. 「対応なし・標本数が異なる」場合の多重比較


    分散分析で有意差がなければ,下位検定としての水準間を比較する多重比較は行わない.ここでは練習のために.


    もしも有意差があれば,どの水準の触媒が安く製造できるか知りたくなる.

    しかし分散分析は,どの水準とどの水準に差があるかは,明らかにできない.

    このため多重比較という事後検定(下位検定)に進み,水準(群)間の比較を行う.

    Tukey(テューキー)の方法は各水準の標本数が数が同じ場合に用いる.

    標本数が異なる場合はテューキー・クレーマーの方法 Tukey-Kramer's method を用いる.



    1. 「対応なし・標本数が異なる」場合の青木の Tukey-Kramer(テューキー・クレーマー)の方法

      Rの「コンソール」画面に,

      source("http://aoki2.si.gunma-u.ac.jp/R/src/tukey.R", encoding="euc-jp") # を貼り付け青木先生のプログラムとリンクさせる.

      # 心配性の受講生への青木先生のコメント

      tukey ( y$製造原価, y$触媒 ) # と記述し,青木のTukey-Kramer(テューキー・クレーマー)の方法で多重比較を行う.

      > source("http://aoki2.si.gunma-u.ac.jp/R/src/tukey.R", encoding="euc-jp")
      > tukey ( y$製造原価, y$触媒 )
      $result1
             n     Mean Variance
      Group1 8 16676.25 405960.5
      Group2 4 15892.25 117901.6
      Group3 3 16609.33 131656.3
      
      $Tukey
                  t          p
      1:2 2.3846906 0.08172134
      1:3 0.1841091 0.98151564
      2:3 1.7488108 0.22797705
      
      $phi
      [1] 12
      
      $v
      [1] 288228.4
      


      分析結果

      多重比較の帰無仮説は二群間の平均値には差が無いである.

      [$Tukey]の欄で有意水準5%以下のp値を探すと,どの水準間にも有意差のある触媒はなかった.



    2. 「対応なし・標本数が異なる」場合の テューキーの方法(Tukey's honestly significant difference test)

      標本数が異なる場合はテューキー・クレーマーの方法 Tukey-Kramer's method を用いる.

      しかし,Rの関数「TukeyHSD」(テューキーの方法) は標本数が異なった場合も使えるようにTukey-Kramerのへ拡張が行われているかもしれない?.

      > help(TukeyHSD)

      "Technically the intervals constructed in this way would only apply to balanced designs where there are the same number of observations made at each level of the factor.

      This function incorporates an adjustment for sample size that produces sensible intervals for mildly unbalanced designs."



      Rの「コンソール」画面に,

      TukeyHSD ( aov ( y$製造原価 ~ y$触媒) )   と記述し,Tukey(テューキー)の方法で多重比較を行う.

      require(graphics) # グラフパッケージの読み込み

      plot (TukeyHSD ( aov ( y$製造原価 ~ y$触媒) ) ) # 結果をグラフで表す


      > TukeyHSD ( aov ( y$製造原価 ~ y$触媒) ) 
        Tukey multiple comparisons of means
          95% family-wise confidence level
      
      Fit: aov(formula = y$製造原価 ~ y$製造原価)
      
      $`y$触媒`
                          diff        lwr        upr       p adj
      水準A2-水準A1  -784.00000  -1661.0971    93.09706  0.0817213
      水準A3-水準A1   -66.91667  -1036.5840   902.75062  0.9815156
      水準A3-水準A2   717.08333   -376.8489  1811.01557  0.2279770
      



      分析結果

      多重比較の帰無仮説は二群間の平均値には差が無いである.

      [p adj]の欄で有意水準5%以下を探すと,どの水準間にも有意差のある触媒はなかった.

      グラフの見方は横棒グラフが,ゼロの値の位置をまたがなければ,有意であると読む.


    3. 「対応なし・標本数が異なる」場合の ボンフェローニ(Bonferroni)法

      Rの「コンソール」画面に,

      pairwise.t.test( y$製造原価, y$触媒, p.adjust.method="bonferroni" )

      と記述し,Bonferroni法で多重比較を行う.

      >  pairwise.t.test( y$製造原価, y$触媒, p.adjust.method="bonferroni" )
      
              Pairwise comparisons using t tests with pooled SD 
      
      data:  y$製造原価 and y$触媒 
      
             水準A1 水準A2
      水準A2 0.10   -     
      水準A3 1.00   0.32  
      
      P value adjustment method: bonferroni

      分析結果

      多重比較の帰無仮説は二群間の平均値には差が無いである.有意水準5%以下を捜すと,どの水準間にも,有意な差がないことがわかる.


    4. 「対応なし・標本数が異なる」場合の Holm法

      Rの「コンソール」画面に,

      pairwise.t.test( y$製造原価, y$触媒, p.adjust.method="holm" )

      と記述し,Holm法で多重比較を行う.

      >  pairwise.t.test( y$製造原価, y$触媒, p.adjust.method="holm" )
      
              Pairwise comparisons using t tests with pooled SD 
      
      data:  y$製造原価 and y$触媒 
      
             水準A1 水準A2
      水準A2 0.10   -     
      水準A3 0.86   0.21  
      
      P value adjustment method: holm 
      

      分析結果

      多重比較の帰無仮説は二群間の平均値には差が無いである.有意水準5%以下を捜すと,どの水準間にも,有意な差がないことがわかる.

一元配置分散分析と多重比較 (対応なし・標本数が異なる)     [ 目次に戻る ]