タグ:R ( 3 ) タグの人気記事

記事の下の方に広告が表示される場合があります。この広告はエキサイトの広告枠です。
 

R 実験計画法 直交表の分散分析(excelデータの活用)

はじめに

excelで直交表を作成し,クリップボードを経由してRで分散分析をする手順です。私の環境がWindows XPなので,別のOSの場合は違うかもしれません。
 扱うデータはL9の直交表です。実際の実験では割り当てられる因子の数も少なく、交互作用の影響が大きくでた時には問題もあるので、よく分かっている方を除いてはL9はお勧めしません。解説用にサイズが小さいので良く利用されることがありますが、あまりお勧めしません。


関数の定義

 この作業を始める前に以下の関数を定義して下さい。関数の名前は適当に変えてくださいね。

関数の定義


##excelのデータをクリップボード経由で取込む
excel <- function(){
read.delim("clipboard")
}

## 直交表を因子に変換
## 最後の列にデータとすること
factor2 <-function(data){
  n=ncol(data)-1
  for(i in 1:n){
    data[,i]<-factor(data[,i])
  }
  return(data)
}




手順

 直交表を扱うときに,データフレームの水準を数値ではなく因子に変換する作業が必要になります。今回定義したfactor2は指定データフレームの指定列を数値から因子に変換します。この作業で最も重要となるのがこの数値から因子への変換です。

1. excelで直交表を作成し,データを入力する。
 この時に、実験データは最後の列に記入して下さい。関数定義で最後の列以外を因子に変換するようにしています。

a0031863_21463471.jpg


2. 直交表の列名を含めてクリップボードに入れる。(コピーをする)
a0031863_21465250.jpg


3. 以下のコマンドを入力し,dataにクリップボードのデータを入れる。

コマンド


data <- excel()



4.以下のコマンドを入力し,直交表の水準データを数値から因子にする。

コマンド


data <- factor2(data)
c(data)



a0031863_226538.jpg


このときに,因子は Levels: 1 2 3 などと表示されますのでちゃんと因子に変換できていることを確認して下さいね。間違っていると分散分析結果がおかしくなります。

5. 分散分析をします。

分散分析


summary(aov(y~A+B+C,data=data))



a0031863_21472352.jpg

※ちゃんと自由度が計算されているのかチェックして下さいね。自由度が間違っている場合は数値から因子に変換できていません。


最後に

なんだか少しづつ、分かってきた気がします。まだ私は始めたばかりのひよこなので、いろいろと至らないこともあるかと思いますが、少しづつ勉強していきたいと思っています。プログラムについてもまだよく分かっていないので、もう少し効率の良い方法もあるかもしれませんね。


[PR]
記事の下の方に広告が表示される場合があります。この広告はエキサイトの広告枠です。
by potto37 | 2007-12-21 21:45 | Trackback | Comments(0)
 

R 実験計画法 直交表から要因効果の計算

はじめに

統計処理環境であるRを使い始めたら、すごく便利なことに気がつきました。今回は前回の続きのである、直交表から要因効果を計算してみたいと思います。思いっきり、個人メモです。質問に答えられないので予めご了承下さいね。


データについて

データは前回の分散分析に使ったデータです。
ぽっとの陽だまり研究室 : R 実験計画法 直交表の分散分析メモ


データ


data <- data.frame(
A = c("1","1","1","2","2","2","3","3","3"),
B = c("1","2","3","1","2","3","1","2","3"),
C = c("1","2","3","2","3","1","3","1","2"),
e = c("1","2","3","3","1","2","2","3","1"),
y = c(1.0,1.2,1.3,1.1,1.0,1.4,1.4,1.5,1.6)
)



簡単な解説

前回のL9のデータから各水準の効果について調べたいと思います。まず、直交表の各因子について、各水準毎の平均を求めます。その求めた平均が各因子の水準の効果です。

例えば、因子Aについて水準1の効果は、因子Aの水準1の平均になるので、以下のコマンドで因子Aの水準1の効果を得ることが出来ます。水準2,3についても同様にコマンドを与えることで得ることが出来ます。

例1



## 因子Aの水準1の効果
mean(subset(data,A==1)$y)

## 因子Aの水準2の効果
mean(subset(data,A==2)$y)

## 因子Aの水準3の効果
mean(subset(data,A==3)$y)




簡単なプログラムで

直交表が大きくなると、上記方法では面倒なので簡単なプログラムを作ってみます。dataからdatに列が因子、行が水準で計算結果を代入しています。もっと良い方法があると思いますが、とりあえず簡単な方法で書いてみました。

簡単なプログラム例


## 要因効果の抽出
dat <- matrix(, nrow=3, ncol=3)
for(j in 1:3){
  for(i in 1:3){
   dat[i,j] = mean(subset(data,data[,j]==i)$y)
  }
}
dat



最後に

Rって便利ですね。以前は、重回帰分析のプログラムを自分で作ってデータの分析をしたこともあったけど、Rならば簡単に重回帰分析もできて、いろいろな統計処理が簡単に出来るので便利です。そしてなによりも、多くの研究者達が利用されているので、情報を得ることができるのも良いですね。


[PR]
記事の下の方に広告が表示される場合があります。この広告はエキサイトの広告枠です。
by potto37 | 2007-12-16 18:30 | Trackback | Comments(0)
 

R 実験計画法 直交表の分散分析メモ

はじめに

まず始めに、私はRについては今始めたばかりの初心者です。普段は手計算で実験データの解析をしていたのですが、さすがに面倒になったので、統計処理の環境であるRを始めてみました。Rの本を注文したのですが、何故か1週間経ってもこないので、Webで勉強してとりあえず計算しています。
 このメモは個人的なメモで質問には答えられるようなレベルでないことを予めご了承下さい。



Rで分散分析(L9直交表)


 簡単な例題として、L9の直交表についての計算例です。実際の実験ではL9ではいろいろと問題もあるかと思いますが、例題としてはサイズが小さいので選びました。データは書籍の「品質管理のための 実験計画法テキスト」から利用させて頂きました。(書籍と計算結果を比較すると楽なので利用させて頂きました。)


直交表の分散分析例


## 中里博明ほか「品質管理のための 実験計画法テキスト」232頁(日科技連 1997年)

## 直交表
data <- data.frame(
A = c("1","1","1","2","2","2","3","3","3"),
B = c("1","2","3","1","2","3","1","2","3"),
C = c("1","2","3","2","3","1","3","1","2"),
e = c("1","2","3","3","1","2","2","3","1"),
y = c(1.0,1.2,1.3,1.1,1.0,1.4,1.4,1.5,1.6)
)

## 直交表のデータを表示
data

## 分散分析
summary(aov(y ~ A+B+C,data=data))


## 分散分析(Cをeにプーリング)
summary(aov(y ~ A+B,data=data))

## summary(aov(y ~ A+B-C,data=data)) #の方が良いかも




a0031863_22224844.jpg


注意

各因子の水準について文字列として設定する必要があります。(このことに気がつかなくて、苦労しました。)
excelからクリップボードを利用してデータを取り込む場合も、a1,a2,a3などと文字列で直交表を作るとちゃんと計算できるみたいです。
ちゃんと計算できていないと、自由度が違うのですぐに気がつくと思いますけど・・・


追記


数値を因子に変換するfactorコマンドを用いると、上記問題を解消できます。詳しくは以下のリンクをご覧下さいね。

ぽっとの陽だまり研究室 : R 実験計画法 直交表の分散分析(excelデータの活用)



最後に

Rで統計処理はじめたら、意外と便利であることに気がつきました。ちょっとしたプログラムを作ることもできるので、いろいろな処理ができて、すごく便利ですね。手計算でやるよりも速いので嬉しいです。(当たり前だよね)


[PR]
記事の下の方に広告が表示される場合があります。この広告はエキサイトの広告枠です。
by potto37 | 2007-12-11 22:23 | Trackback | Comments(0)
掲載の記事や画像などすべての複写・転載・公衆送信等はご遠慮ください。