Last updated: January 2010. Kajiyama                    [ 目次に戻る ]

原料4種類,反応温度5種類での収率分析: 二元配置の分散分析と多重比較

(因子Aに対応なしX因子Bに対応なし・標本数は1)(多重比較・交互作用なし)


  1. 「標本数は1」の場合の課題

    「繰り返しのない,標本数1」の二元配置のデータを分散分析することは,

    総変動から,因子A,因子Bによる変動と誤差変動を分離して,測定値への影響を検討することである.

    繰り返しの実験を行っていないので,因子間の交互作用AXBは検討できない.

    繰り返しのない実験を行えるのは,因子Aと因子B間に交互作用が存在しないとわかっている場合だけに行う実験計画である.



    作業を改善するために実験を行った.その結果表1の収率の20回の実験データを得た.

    実験は,原料(因子A 4水準)と反応温度(因子B 5水準)の各水準を組み合わせた4x5=20組の条件に対して,繰り返しをせずに,一回ずつランダムに行った.

    4種類の原料間で収率に差があるか,また,5種類の反応温度間で収率に差があるか検定する.



    表1  原料4種類,反応温度5種類での収率%

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

    収率(yield)は生産工程で、理論上可能な物質の最大量(理論収量)に対する実際に得た量(収量)の比率である.

    収率はその工程がすぐれているかどうかの指標の一つである.


     
    因子A 原料
    水準A1 水準A2 水準A3 水準A4


    B



    水準B1 080゜C 89.4 89.7 89.3 89.9
    水準B2 090゜C 90.3 90.0 90.1 89.7
    水準B3 100゜C 91.8 91.5 90.8 90.1
    水準B4 110゜C 92.1 91.9 91.6 92.3
    水準B5 120゜C 93.4 93.3 92.8 92.9

    測定回数は1回である.測定のくり返しはない.


    表2  原料4種類,反応温度5種類での収率%

    原料
    (因子A)
    A1 A2 A3 A4
    反応温度
    (因子B)
    080
    ゜C
    090
    ゜C
    100
    ゜C
    110
    ゜C
    120
    ゜C
    080
    ゜C
    090
    ゜C
    100
    ゜C
    110
    ゜C
    120
    ゜C
    080
    ゜C
    090
    ゜C
    100
    ゜C
    110
    ゜C
    120
    ゜C
    080
    ゜C
    090
    ゜C
    100
    ゜C
    110
    ゜C
    120
    ゜C
    収率%
    (測定値)
    89.4 89.7 89.3 89.9 90.3 90.0 90.1 89.7 91.8 91.5 90.8 90.1 92.1 91.9 91.6 92.3 93.4 93.3 92.8 92.9

    測定回数は1回である.測定のくり返しはない.


    表3 原料4種類,反応温度5種類での収率%

    因子B 因子A 測定値
    反応温度 原料 収率
    080゜C A1 89.4
    080゜C A2 89.7
    080゜C A3 89.3
    080゜C A4 89.9
    090゜C A1 90.3
    090゜C A2 90.0
    090゜C A3 90.1
    090゜C A4 89.7
    100゜C A1 91.8
    100゜C A2 91.5
    100゜C A3 90.8
    100゜C A4 90.1
    110゜C A1 92.1
    110゜C A2 91.9
    110゜C A3 91.6
    110゜C A4 92.3
    120゜C A1 93.4
    120゜C A2 93.3
    120゜C A3 92.8
    120゜C A4 92.9

    変数の「反応温度と原料」因子のデータ形式はfactorである.

        表4 原料4種類,反応温度5種類での収率% 「スタック形式」

    反応温度 原料 収率
    080゜C A1 89.4
    080゜C A2 89.7
    080゜C A3 89.3
    080゜C A4 89.9
    090゜C A1 90.3
    090゜C A2 90.0
    090゜C A3 90.1
    090゜C A4 89.7
    100゜C A1 91.8
    100゜C A2 91.5
    100゜C A3 90.8
    100゜C A4 90.1
    110゜C A1 92.1
    110゜C A2 91.9
    110゜C A3 91.6
    110゜C A4 92.3
    120゜C A1 93.4
    120゜C A2 93.3
    120゜C A3 92.8
    120゜C A4 92.9

    変数の「反応温度と原料」因子のデータ形式は factorである.



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


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

    > x <- read.table ("clipboard", header=TRUE )
    > x 
       反応温度 原料 収率
    1    080゜C   A1 89.4
    2    080゜C   A2 89.7
    3    080゜C   A3 89.3
    4    080゜C   A4 89.9
       ・・・・・・
    16   110゜C   A4 92.3
    17   120゜C   A1 93.4
    18   120゜C   A2 93.3
    19   120゜C   A3 92.8
    20   120゜C   A4 92.9
    
    > str (x) 
    > str (x)
    'data.frame':   20 obs. of  3 variables:
     $ 反応温度: Factor w/ 5 levels "080゜C","090゜C",..: 1 1 1 1 2 2 2 2 3 3 ...
     $ 原料    : Factor w/ 4 levels "A1","A2","A3",..: 1 2 3 4 1 2 3 4 1 2 ...
     $ 収率    : num  89.4 89.7 89.3 89.9 90.3 90 90.1 89.7 91.8 91.5 ...
    


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

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

    表4 原料4種類,反応温度5種類での収率%「スタック形式」の緑の部分をマウスで選択し,これをコピーする.

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

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

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

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

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


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

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


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

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

    Rが理解した,サンプル・分析値・水準の3変数名のリストと変数ごとのカテゴリの数とカテゴリ名が示される.

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


    注意: 因子の反応温度と原料のデータ形式が [Factor] 形式でないときは,

    x$反応温度 <- factor( x$反応温度 )

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


  3. 1因子ごとの基礎統計


    a. 因子AB別の基礎統計量

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

    by ( x$収率, x$反応温度, summary)

    by ( x$収率, x$原料, summary) #   と記述し,「スタック形式」のデータで,各因子の群別の基礎統計利用を計算する.


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

    plotMeans( x$収率, x$反応温度, error.bars="se")

    plotMeans( x$収率, x$原料, error.bars="se") #  と記述し,原料ごとの平均値の推移図を表示する.


    c. グラフの観察




  4. 2因子間の基礎統計 2因子間AXBの交互作用の可能性検討


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

    どの因子間の水準の組み合わせの時に,収率が高いかがわかる.ただし標本数が1の場合は検定できない.

    plotMeans( x$収率, x$反応温度, x$原料, error.bars="se")


    b. グラフの観察




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


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

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

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


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

    bartlett.test ( x$収率 ~ x$原料 )

    bartlett.test ( x$収率 ~ x$反応温度 )

    と記述し,「スタック形式」のデータで,分散の等質性の検定を行う.

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

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


    >  bartlett.test ( x$収率 ~ x$原料 )
    
            Bartlett test of homogeneity of variances
    
    data:  x$収率 by x$原料 
    Bartlett's K-squared = 0.0839, df = 3, p-value = 0.9937
    
    > bartlett.test ( x$収率 ~ x$反応温度 )
    
            Bartlett test of homogeneity of variances
    
    data:  x$収率 by x$反応温度 
    Bartlett's K-squared = 5.6648, df = 4, p-value = 0.2256
    

    分析結果

    p値が0.99と0.23であり,p値>0.05であるので2つの因子の各水準の分散は等しいと結論する.



  6. 「因子Aに対応なしX因子Bに対応なし・標本数は1」aov関数による二元配置分散分析

    「因子Aに対応あり・標本数が同じ」1要因の分散分析と比べてみる



    4種類の原料間で収率に差があるか,また,5種類の反応温度間で収率に差があるか,有意水準5%で検討する.

    帰無仮説は「各水準の母集団の平均値は等しい」.

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



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

    summary ( aov ( 収率 ~ 反応温度+原料, data=x ) )

    と記述し,「スタック形式」のデータで,「繰り返しのない二元配置の分散分析」を行う.

    >  summary ( aov ( 収率 ~ 反応温度+原料, data=x ) )
                Df Sum Sq Mean Sq F value    Pr(>F)    
    反応温度     4 32.957   8.239 52.9572 1.586e-07 ***
    原料         3  0.806   0.269  1.7258    0.2148    
    Residuals   12  1.867   0.156                      
    ---
    Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
    
    
    分散分析表が表示される.これを整形すると,以下のようになる.
    -----------------------------------------------------------          (変動) (不偏分散) (分散比) 変動要因 自由度 平方和 平均平方和 F値 P値 ----------------------------------------------------------- 反応温度 4 32.957 8.239 52.9572 0.000 *** 原料 3 0.806 0.269 1.7258 0.215 残差 12 1.867 b0.156 ----------------------------------------------------------- *** p<.001
    分析結果

    有意水準5%で検定を行ったところ,反応温度の差においてはP値<.05であるから,帰無仮説は棄却する.

    すなわち反応温度による収率の差はあるとする.

    原料の差においてはP値>.05であるから帰無仮説を採択する.すなわち原料による収率の差があるとはいえない.

    反応温度は,収率に影響を与えている.



  7. 「交互作用なし・因子Aに対応なしX因子Bに対応なし・標本数は1」 青木先生のtwoway.anova関数による分散分析


    4種類の原料間で収率に差があるか,また,5種類の反応温度間で収率に差があるか,有意水準5%で検討する.

    帰無仮説は「各水準の母集団の平均値は等しい」.

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

    注意:青木のtwoway.anova関数では,因子の水準に1,2,3,4,5,・・・の数字を用いる.表3の水準を数字に書き換えたのが表4である


    表4 原料4種類,反応温度5種類での収率% 「スタック形式」

    反応温度 原料 収率
    1 1 89.4
    1 2 89.7
    1 3 89.3
    1 4 89.9
    2 1 90.3
    2 2 90.0
    2 3 90.1
    2 4 89.7
    3 1 91.8
    3 2 91.5
    3 3 90.8
    3 4 90.1
    4 1 92.1
    4 2 91.9
    4 3 91.6
    4 4 92.3
    5 1 93.4
    5 2 93.3
    5 3 92.8
    5 4 92.9


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

    > y <- read.table ("clipboard", header=TRUE)
    >
    >  str (y)
    'data.frame':   20 obs. of  3 variables:
     $ 反応温度: int  1 1 1 1 2 2 2 2 3 3 ...
     $ 原料    : int  1 2 3 4 1 2 3 4 1 2 ...
     $ 収率    : num  89.4 89.7 89.3 89.9 90.3 90 90.1 89.7 91.8 91.5 ...
    
    >  y$反応温度 <- factor( y$反応温度 )
    >  y$原料 <- factor( y$原料 ))
    >  str (y)
    'data.frame':   20 obs. of  3 variables:
     $ 反応温度: Factor w/ 5 levels "1","2","3","4",..: 1 1 1 1 2 2 2 2 3 3 ...
     $ 原料    : Factor w/ 4 levels "1","2","3","4": 1 2 3 4 1 2 3 4 1 2 ...
     $ 収率    : num  89.4 89.7 89.3 89.9 90.3 90 90.1 89.7 91.8 91.5 ...
    

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

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

    表4 原料4種類,反応温度5種類での収率%「スタック形式」の緑の部分をマウスで選択し,これをコピーする.

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

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

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

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

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

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

    Rが理解した,サンプル・分析値・水準の3変数名のリストと変数ごとのカテゴリの数とカテゴリ名が示される.

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


    注意: 因子の反応温度と原料のデータ形式が [Factor] 形式でないときは,

    y$反応温度 <- factor( y$反応温度 )

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

    (4) 分散分析

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

    source("http://aoki2.si.gunma-u.ac.jp/R/src/twoway_anova.R", encoding="euc-jp")

    twoway.anova(y$収率, y$原料, y$反応温度)

    と記述し,「スタック形式」のデータで,二元配置の分散分析を行う.

    > source("http://aoki2.si.gunma-u.ac.jp/R/src/twoway_anova.R", encoding="euc-jp")
    > twoway.anova(y$収率, y$原料, y$反応温度)
           SS d.f.        MS   F value      P value
    A  0.8055    3 0.2685000  1.725763 2.147822e-01
    B 32.9570    4 8.2392500 52.957151 1.585848e-07
    e  1.8670   12 0.1555833        NA           NA
    T 35.6295   19 1.8752368        NA           NA
    
    
    分散分析表が表示される.これを整形すると,以下のようになる.
    -----------------------------------------------------------          (変動) (不偏分散) (分散比) 変動要因 自由度 平方和 平均平方和 F値 P値 ----------------------------------------------------------- 原料 3 0.806 0.269 1.726 0.215 反応温度 4 32.957 8.239 52.957 0.000 *** 残差 12 1.867 0.156 ----------------------------------------------------------- *** p<.001
    分析結果

    有意水準5%で検定を行ったところ,反応温度の差においてはP値<.001であるから,帰無仮説は棄却する.

    すなわち反応温度による収率の差はあるとする.

    原料の差においてはP値>.05であるから帰無仮説を採択する.すなわち原料による収率の差があるとはいえない.

    反応温度は,収率に影響を与えている.



  8. 「交互作用なし・因子Aに対応なしX因子Bに対応なし・標本数は1」の二元配置の多重比較

     
    二元配置で交互作用がなかったので,因子A・因子Bについて多重比較を検討すればよい.交互作用があれば,因子A・因子B の多重比較は不要になる.

    因子B「反応温度」の水準間に差があることが分散分析によりわかった.

    分散分析では,「反応温度」のどの水準とどの水準に差があるかは,明らかにしてくれない.

    このため多重比較の事後検定(下位検定)を行う.


    1. 「対応なしX対応なし・標本数は1」場合の テューキーHSD(Tukey's‘Honest Significant Difference’method)の方法では

      関数TukeyHSDの使用法は,Rの「コンソール」画面に,

      ? TukeyHSD

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

      TukeyHSD(aov ( 収率 ~ 反応温度 + 原料, data=x ), "反応温度", ordered = TRUE) # 分散分析の結果をTukeyHSD( )に入れる

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

      plot(TukeyHSD(aov ( 収率 ~ 反応温度 + 原料, data=x ), "反応温度")) # TukeyHSDの結果をグラフにする

      plot(TukeyHSD(aov ( 収率 ~ 反応温度 + 原料, data=x ), "原料")) # TukeyHSDの結果をグラフにする

      と記述し,「データセットx」についてTukeyHSD(テューキー)の方法で多重比較を行う.

      > TukeyHSD(aov ( 収率 ~ 反応温度 + 原料, data=x ), "反応温度", ordered = TRUE) #「データセットx」
      
        Tukey multiple comparisons of means
          95% family-wise confidence level
          factor levels have been ordered
      
      Fit: aov(formula = 収率 ~ 反応温度 + 原料, data = x)
      
      $反応温度
                     diff         lwr      upr     p adj
      090゜C-080゜C 0.450 -0.43901175 1.339012 0.5166792
      100゜C-080゜C 1.475  0.58598825 2.364012 0.0014637
      110゜C-080゜C 2.400  1.51098825 3.289012 0.0000145
      120゜C-080゜C 3.525  2.63598825 4.414012 0.0000002
      100゜C-090゜C 1.025  0.13598825 1.914012 0.0216926
      110゜C-090゜C 1.950  1.06098825 2.839012 0.0001160
      120゜C-090゜C 3.075  2.18598825 3.964012 0.0000010
      110゜C-100゜C 0.925  0.03598825 1.814012 0.0401310
      120゜C-100゜C 2.050  1.16098825 2.939012 0.0000712
      120゜C-110゜C 1.125  0.23598825 2.014012 0.0117273
      
      > require(graphics) # グラフパッケージの読み込み
      > plot(TukeyHSD(aov ( 収率 ~ 反応温度 + 原料, data=x ), "反応温度")) # TukeyHSDの結果をグラフにする
      > plot(TukeyHSD(aov ( 収率 ~ 反応温度 + 原料, data=x ), "原料"))



      分析結果

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

      グラフでは95%の信頼区間が0を含まなければ有意な差があることを示す.分散分析結果が示したように原料間は全て有意ではない.

      [p adj]の欄で有意水準5%以下を捜すと,[090゜C-080゜C p=0.5166792]を除く組み合わせで有意差がある.

      120゜Cのとき収率は最大で,090℃〜080℃のとき最小になる.


    2. 「対応なしX対応なし・標本数は1」の場合の ボンフェローニ(Bonferroni)法

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

      pairwise.t.test(x$収率 , x$反応温度, p.adjust.method="bonferroni" )

      pairwise.t.test(x$収率 , x$原料, p.adjust.method="bonferroni" )


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

      > pairwise.t.test(x$収率 , x$反応温度, p.adjust.method="bonferroni" )
      
              Pairwise comparisons using t tests with pooled SD 
      
      data:  x$収率 and x$反応温度 
      
             080゜C  090゜C  100゜C  110゜C
      090゜C 1.0000  -       -       -     
      100゜C 0.0018  0.0369  -       -     
      110゜C 8.1e-06 9.5e-05 0.0733  -     
      120゜C 5.4e-08 3.4e-07 5.3e-05 0.0186
      
      P value adjustment method: bonferroni 
      >
      >
      > pairwise.t.test(x$収率 , x$原料, p.adjust.method="bonferroni" )
      
              Pairwise comparisons using t tests with pooled SD 
      
      data:  x$収率 and x$原料 
      
         A1 A2 A3
      A2 1  -  - 
      A3 1  1  - 
      A4 1  1  1 
      
      P value adjustment method: bonferroni 
      

      分析結果

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

      水準B1〜B5間で有意水準5%以下を探すと,

      [080℃-090]と[100℃-110]を除いた組み合わせに有意差がある.

      この結果,反応温度120℃で収率は最大になり,反応温度090℃〜080℃で収率は最小になる.


    3. 「対応なしX対応なし・標本数は1」の場合の Holm法

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

      pairwise.t.test(x$収率 , x$反応温度, p.adjust.method="holm" )

      pairwise.t.test(x$収率 , x$原料, p.adjust.method="holm" )


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

      > pairwise.t.test(x$収率 , x$反応温度, p.adjust.method="holm" )
      
              Pairwise comparisons using t tests with pooled SD 
      
      data:  x$収率 and x$反応温度 
      
             080゜C  090゜C  100゜C  110゜C 
      090゜C 0.15240 -       -       -      
      100゜C 0.00089 0.01107 -       -      
      110゜C 6.5e-06 5.7e-05 0.01466 -      
      120゜C 5.4e-08 3.0e-07 3.7e-05 0.00743
      
      P value adjustment method: holm 
      >
      >
      > pairwise.t.test(x$収率 , x$原料, p.adjust.method="holm" )
      
              Pairwise comparisons using t tests with pooled SD 
      
      data:  x$収率 and x$原料 
      
         A1 A2 A3
      A2 1  -  - 
      A3 1  1  - 
      A4 1  1  1 
      
      P value adjustment method: holm 
      

      分析結果

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

      水準B1〜B5間で有意水準5%以下を探すと,

      [080℃-090]を除いた組み合わせに有意差がある.

      この結果,反応温度120℃で収率は最大になり,反応温度090℃〜080℃で収率は最小になる.

くり返しのない二元配置分散分析 (標本数は1)     [ 目次に戻る ]