os xでキーボードのリピートを爆速にする

ubuntuやwindowsからmacに移行してから、keyboardのリピートが遅くていらだったことはありませんか？ keyremap4macbookというソフトウェアを用いれば、その悩みは解消されます。インストールが完了したらアプリを起動して、Key Repeatのタブをクリックして、お好きなように、リピート速度を決定してください。

https://pqrs.org/macosx/keyremap4macbook/index.html.ja

VirtualboxのゲストのUbuntuにポートフォワードで接続する

以前、virtualboxにインストールしたUbuntuに手の込んだ方法で接続する方法を紹介しました

今回は、ポートフォワードという方法で接続する紹介を勉強しましたので、記事にします。

設定したいOSのSettingsをクリックしてNetworkの設定のタブに移動します。そして、Adapter 1の設定画面でPort Forwardingをクリックします。

すると以下のような画面が出ます。右はしのプラスボタンをクリックして以下のようにSSHとhttpのポートを開放します。

以上で設定完了です。

SSH接続は以下のように行います。

$ ssh kappa@localhost -p 2222

ゲストのUbuntuを起動してapacheが起動した状態で、以下のようにしてhttp://localhost8080に接続すれば、サーバー側のディレクトリを閲覧することが可能です。

以下のサイトを参考にさせて頂きました。
http://www.karakaram.com/virtualbox-port-fowarding

VirtualboxのゲストOSの環境をコピーする

Virtualboxを使用していると、ついついチャレンジングなことをゲスト環境でしてしまい、復旧に手間がかかることがあります。そういった時にゲストOSの環境をそのままコピーしておくことができると大変便利です。それには以下のようにターミナルを用いた処理を必要とします。以下、VirtualboxにUbuntuというゲストOSが入っている前提で話を進めます。

#以下のディレクトリに移動する

$ pwd

/Users/kappa/VirtualBox VMs/Ubuntu

$ ls

Logs             Ubuntu.vbox      Ubuntu.vdi

Snapshots        Ubuntu.vbox-prev 

#virtualboxに備え付けのコピーコマンドを実行する

$ VBoxManage clonehd Ubuntu.vdi Ubuntu_copy.vdi

$ ls

Logs             Ubuntu.vbox      Ubuntu.vdi

Snapshots        Ubuntu.vbox-prev Ubuntu_copy.vdi

これでゲスト環境のコピーファイルが完成しました。Ubuntu_copy.vdiです。最後にこれをVirtualboxに読み込みます。以下はGUI作業です。

最初はこんな感じです。

 FileからVirtual Media Managerを開きます。

 Copyをクリックします。

先ほどターミナルでコピーしたUbuntu_copy.vdiを選択して後はデフォルトの設定に任せていけば、仮想環境がコピーされます。ただし、この段階では仮想ハードディスクがコピーされただけみたいです。

最後に普段行うように、新しいゲストOSの作成を行います。

ここがミソです。Use an existing virtual hard drive fileを選択します。そして、先ほどこコピーしたものを選択します。

これでコピー終了です。

ちゃんと起動し、ログイン画面が表示されました。

Virtualbox4.3.2にWindows8.1をインストールする

Virtualbox4.3.2にWindows8.1(64bit)をインストールをしたときのメモを書きます。

ポイントは2つです。
一つ目はOSが64bitなら、virtualbox側も64bitにそろえること。
二つ目はVirtualbox4.3.2にWindows8.1(64bit)の組み合わせでは、インストール時にエラーが出てしまいます。そこで以下のサイトを参考にして、仮想ハードディスク（名前はWindows8とした）ターミナルからvirtualboxに向けてコマンドを入力してやる必要があります。
$ VBoxManage setextradata Windows8 VBoxInternal/CPUM/CMPXCHG16B 1

参考サイト
http://zak-za-k.blogspot.jp/2013/07/virtualbox-for-macwindows81error-code.html
http://fnya.cocolog-nifty.com/blog/2013/06/virtualbox-wi-1.html

作業時のスクリーンキャプチャを添付しておきます。
以下のようなバージョンで仮想ハードディスクを設定します。

仮想ハードディスクを作成した後はこんな感じ。

killを使ってアプリケーションを強制終了

OS Xを使用していると、アプリケーションの強制終了に苦慮することがあります。そこで今回は、killコマンドを用いて指定したアプリケーションを強制終了することを考えます。

今、Google Chromeがフリーズしてしまったとします。
こういった時は、ターミナルを起動して、以下の作業を行います。

ちなみにpsのカラムは、以下のようになります。

PID : Process IDentification
TTY : そのプロセスが実行されているターミナルの生。??であれば、そのプロセスがGUI画面以外で実行されていことを意味する。
TIME : 特定のプロセスを実行するのに要した時間
CMD : 実行しているプロセスの完全なパス名

【参考文献】
Dave Taylor著　酒井皇治訳『入門 Unix for OS X 』オライリージャパン 2013

LaTeXItの紹介

Texを使用すると数式がとても奇麗に描写できます。ただ、ソースコードを書いてコンパイルするのは若干煩雑です。特に、たった２−３行の数式が必要なだけなのに、バグに悩まされるのはちょっとつらい時があります。そこで、今回はLaTeXItというアプリケーションを紹介します。
インストールはhttp://tug.org/mactex/のサイトからMacTeXをインストールするだけです。自動的に使用可能となります。

上のウィンドウの下半分にあるスペースにtexのコマンドを入力して、右下のLaTeX it!をクリックすれば上半分に数式が表示されます。数式の上で右クリックをすると、様々な形式で数式をコピーすることができます。

command + l というショートカットキーで右下のLaTeX it!をクリックを省略することを覚えておくと便利です。

病理のバーチャルスライド

最近はインターネットが発達し、手術の摘出標本から作成したスライドをネットで見ることが可能となってきました。

一例を紹介します。

Virtual Slide Box of  Histology

http://www.path.uiowa.edu/virtualslidebox/

http://www.path.uiowa.edu/virtualslidebox/nlm_histology/inde

逐一、サーバーとデータ通信を行うため、ネット環境によっては動作が非常に緩慢ですが、教科書よりも間違いなく臨場感のある勉強ができます。

第一種の誤りとは

第一種の誤りとは

「帰無仮説が真のとき、これを棄却してしまう」

誤りのことです。

そして、第一種の誤りを犯す確率をαとして、有意水準または危険率と一般に呼びます。

第一種の誤りについてもう少し直感的に理解するには、条件付き確率を扱っていることを認識する必要があります。

すなわち、危険率とは

1. 帰無仮説が真のとき（こういう条件の下だと）
2. これを棄却してしまう確率

ということなのです。

ですので、あくまで「帰無仮説が正しいとき」を想定しなければ、危険率αのニュアンスはわからないのです。

正規分布からのサンプルの平均値の検定を題材にして、危険率αとはどのようなものなのかということをシミュレーションを通じて考えてみたいと思います。

今、平均値μ、分散σ^2/μの正規分布(標準正規分布)から、サンプリング(n)を行い、標本平均の分布について考えたいと思います。

一般に平均μ、分散σ^2の正規分布から抽出されたサンプル(n)の平均値の分布は平均値μ、分散σ^2/nの正規分布に従うとされます。

サンプルの平均値を標準化して、

を考えます。Zは標準正規分布に従います。

それでは、平均値12、分散9の正規分布の集団からn=5のサンプリングを行い、その平均値の分布について計算してみます。

一回の試行は以下のようにして行うことができます。

> rnorm(mean=12, sd=3, n=5)

[1] 15.039414  9.264087 11.036560 16.033761 11.038038

Zの値は以下のようにして計算します。

> (mean(rnorm(mean=12, sd=3, n=5) - 12)/(3/sqrt(5)))

[1] 13.86206

ところで、標準正規分布を描写すると以下のようになります。

> curve(dnorm(x), -3, 3)

ここで、上側確率、下側確率を計算してグラフに記入

> qnorm(0.025)

[1] -1.959964

> qnorm(0.975)

[1] 1.959964

> abline(v=qnorm(0.025), col="red")

> abline(v=qnorm(0.975), col="red")

それでは今回想定している平均値12、分散9の正規分布の集団からn=5のサンプリングを繰り返し(300回)行い、その都度Z値を計算しています。これらのZ値のうち、上側下側合わせて5%の区間にどの程度入ってくるか調べます。

curve(dnorm(x), -3, 3)

abline(v=qnorm(0.025), col="green")

abline(v=qnorm(0.975), col="green")

for(i in 1 : 100){

abline(v=(mean(rnorm(mean=12, sd=3, n=5) - 12)/(3/sqrt(5))), col = "magenta")

}

確かに、このように繰り返し平均値をもとめていると、

「帰無仮説は正しい」

のにも関わらず、一定の割合でサンプルの平均値は5%以下の確率で生じるとされる区間に入ってきてしまいます。

従って、ある検定を行いたいと考えたときに、サンプルの統計量が偶然には起こりそうにも無い値であったとしても、帰無仮説が誤っている可能性の他に、危険率αの確率で帰無仮説が正しい場合にも同様の統計量が計算される可能性があるといえるのです。

【参考文献】

山田剛史ほか『Rによるやさしい統計学』オーム社 2008 109-119pp

発生学のおすすめサイト

The University of New South Walesのサイトがとてもまとまっています。
#UNSW Embryology
http://php.med.unsw.edu.au/embryology/index.php?title=Main_Page


発生学を動画で学べます。
#UNSW Embryology - Movie
http://php.med.unsw.edu.au/embryology/index.php?title=Movies

第一種の誤りとは ー危険率α＝5%の意味ー

第一種の誤りとは

「帰無仮説が真のとき、これを棄却してしまう」

誤りのことです。

そして、第一種の誤りを犯す確率をαとして、有意水準または危険率と一般に呼びます。

第一種の誤りについてもう少し直感的に理解するには、条件付き確率を扱っていることを認識する必要があります。

すなわち、危険率とは

1. 帰無仮説が真のとき（こういう条件の下だと）
2. これを棄却してしまう確率

ということなのです。

ですので、あくまで「帰無仮説が正しいとき」を想定しなければ、危険率αのニュアンスはわからないのです。

正規分布からのサンプルの平均値の検定を題材にして、危険率αとはどのようなものなのかということをシミュレーションを通じて考えてみたいと思います。

今、平均値μ、分散σ^2/μの正規分布(標準正規分布)から、サンプリング(n)を行い、標本平均の分布について考えたいと思います。

一般に平均μ、分散σ^2の正規分布から抽出されたサンプル(n)の平均値の分布は平均値μ、分散σ^2/nの正規分布に従うとされます。

サンプルの平均値を標準化して、

を考えます。Zは標準正規分布に従います。

それでは、平均値12、分散9の正規分布の集団からn=5のサンプリングを行い、その平均値の分布について計算してみます。

一回の試行は以下のようにして行うことができます。

> rnorm(mean=12, sd=3, n=5)

[1] 15.039414  9.264087 11.036560 16.033761 11.038038

Zの値は以下のようにして計算します。

> (mean(rnorm(mean=12, sd=3, n=5) - 12)/(3/sqrt(5)))

[1] 13.86206

ところで、標準正規分布を描写すると以下のようになります。

> curve(dnorm(x), -3, 3)

ここで、上側確率、下側確率を計算してグラフに記入

> qnorm(0.025)

[1] -1.959964

> qnorm(0.975)

[1] 1.959964

> abline(v=qnorm(0.025), col="red")

> abline(v=qnorm(0.975), col="red")

それでは今回想定している平均値12、分散9の正規分布の集団からn=5のサンプリングを繰り返し(300回)行い、その都度Z値を計算しています。これらのZ値のうち、上側下側合わせて5%の区間にどの程度入ってくるか調べます。

> curve(dnorm(x), -3, 3)

> abline(v=qnorm(0.025), col="green")

> abline(v=qnorm(0.975), col="green")

> for(i in 1 : 100){

+ abline(v=(mean(rnorm(mean=12, sd=3, n=5) - 12)/(3/sqrt(5))), col = "magenta")

+ }

確かに、このように繰り返し平均値をもとめていると、

「帰無仮説は正しい」

のにも関わらず、一定の割合でサンプルの平均値は5%以下の確率で生じるとされる区間に入ってきてしまいます。

従って、ある検定を行いたいと考えたときに、サンプルの統計量が偶然には起こりそうにも無い値であったとしても、帰無仮説が誤っている可能性の他に、危険率αの確率で帰無仮説が正しい場合にも同様の統計量が計算される可能性があるといえるのです。

【参考文献】

山田剛史ほか『Rによるやさしい統計学』オーム社 2008 109-119pp

統計学的検定の根本的な概念

統計学的検定の根本的な概念は

母集団から抽出したデータから計算した統計量が、次のいずれに属するかを検討することだと思います。

1. 偶然誤差の範囲にある(それくらいのばらつきは生じても不思議ではない)
2. 偶然誤差の範囲を逸脱している(きっと想定した母集団モデルから逸脱させる何か、すなわち要因が働いているはずである)

以上のように考えると、統計学的検定の意味が少しずつ見えてくる気がします。

成書を読むと、1,2と同様の内容が必ず"さらっと"書かれていることに気がつきます。

Rを用いてdot plotを描写

Rを用いてdot plotを描写する方法を勉強しました。
dot plotには平均値と標準偏差を同時に記述することが多いので、今回もその慣例に従いました。

参考にさせていたいたのは下のリンクです。

http://minato.sip21c.org/swtips/MedicalStat-Rev.pdf

まずは、以下のデータをEXELで用意して、tab区切り形式で保存します。

あとはRを開いて以下のコマンドを実行すればよいです。

$ R

#データの読み込み

> dat <- read.table("body_height.txt", header=T)

> dat <- as.data.frame(dat)

> dat

   Male Sex

1   170   1

2   175   1

3   181   1

4   176   1

5   167   1

6   180   1

7   156   2

8   150   2

9   153   2

10  148   2

11  157   2

12  149   2

#Sex列を要因(factor)にする

> dat$Sex <- as.factor(dat$Sex)

> levels(dat$Sex) <- c('M','F')

> dat

   Male Sex

1   170   M

2   175   M

3   181   M

4   176   M

5   167   M

6   180   M

7   156   F

8   150   F

9   153   F

10  148   F

11  157   F

12  149   F

#各列名単独で列の要素にアクセスできるようにする

> attach(dat)

#平均値、標準偏差を計算

> mean <- tapply(Height,Sex,mean)

> sd <- tapply(Height,Sex,sd)

> is <- c(1,2)+0.15

#ドットプロットを描写

> stripchart(dat$Height~Sex,method="jitter",vert=T,ylab="body height(cm)") 

> points(is,mean,pch=18)

> arrows(is,mean-sd,is,mean+sd,code=3,angle=90,length=.1)

勝率と二項検定

サッカーのあるチームAがあるチームに対して、年間を通じて、Bチームと30試合を行い25勝5敗だったとします。

このデータから、チームAはチームBより強いと言えるでしょうか。

言い換えると、チームAのチームBに対する真の勝率は５割以上なのでしょうか。それとも、25勝5敗は誤差の範囲で勝ち越しているだけなのでしょうか。

この問題を解くには、二項検定という手法が用いられます。コインの裏表の出る確率の検定のときにも用いられる手法です。

#勝率が5割だったとした時の、勝利数の確率分布を描写する

> win <- 0:30

> plot(dbinom(win, 30, 0.5), type="h", xlab-"number of win")

#累積分布関数を描写する

> plot(pbinom(win, 30, 0.5), type="h", xlab="number of win")

#95%の水準に赤線を引く

> abline(h=0.95, col="red")

#5%の危険率で二項検定を行う

#方法1

> 1 - pbinom(25, 30, 0.5)

[1] 2.973806e-05

#方法2

> pbinom(25, 30, 0.5, lower.tail=FALSE)

[1] 2.973806e-05

> pbinom(25, 30, 0.5, lower.tail=FALSE) < 0.05

[1] TRUE

> pbinom(25, 30, 0.5, lower.tail=FALSE) < 0.01

[1] TRUE

> binom.test(25,30, 0.5)

 Exact binomial test

data: 25 and 30 

number of successes = 25, number of trials = 30, p-value = 0.0003249

alternative hypothesis: true probability of success is not equal to 0.5 

95 percent confidence interval:

 0.6527883 0.9435783 

sample estimates:

probability of success 

             0.8333333 

以上の計算結果から、チームAのチームBに対する勝率は有意に50%よりも大きいことが言えます。また、その信頼区間は0.65以上0.95以下です。

今回行った計算方法は、帰無仮説であるp=0.5を他の値に変更することによって、二項分布に従う様々な現象(打率、くじ引き etc)を説明することができます。

【参考文献】

山田剛史ほか『Rによるやさしい統計学』オーム社 2008 第12章

重回帰分析におけるモデルの検討

重回帰分析の手法について分子生物学的な事例をもとにして考えてみます。

ある時刻の遺伝子Mの発現量exp_mが、一時間前の遺伝子n及び遺伝子lの発現量(それぞれexp_n, exp_l)によりどのよう説明されるかというモデルを立てます。

具体的には、

というモデルを想定します。

aを切片、b_1、b_2を偏回帰係数、eを残差と呼びます。

それでは、Rを用いて重回帰分析を行ってみます。

#データの入力

> exp_m <- c(45, 46, 47, 57, 78, 89, 99, 89, 101, 108)

> exp_n <- c(102, 106, 109, 120, 124, 130, 146, 178, 190, 196)

> exp_l <- c(56, 57, 58, 50, 52, 61, 43, 63, 64, 60)

#exp_mを目的変数、exp_n, exp_lを説明変数として重回帰分析を行う

> summary(lm(exp_m~exp_n+exp_l))

Call:

lm(formula = exp_m ~ exp_n + exp_l)

Residuals:

    Min 1Q Median 3Q Max 

-11.530 -5.954 -3.900 3.842 24.637 

Coefficients:

                       Estimate   Std. Error     t value    Pr(>|t|)    

(Intercept)        36.3470      34.7881      1.045    0.330842    

exp_n                0.6858         0.1221     5.618     0.000801 ***

exp_l                -1.0022         0.6710    -1.494     0.178918    

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 11.91 on 7 degrees of freedom

Multiple R-squared: 0.822,Adjusted R-squared: 0.7712 

F-statistic: 16.16 on 2 and 7 DF, p-value: 0.002379 

exp_nのp値が有意となっているということは、冒頭の式のb_1が有意に0よりも大きい値であるということができます。また、exp_nを用いれば、各値の全体平均からのばらつきを有意に説明するとも言い換えることができます(ANOVAの考え方と同じです)。一方で、exp_lのp値は有意ではありませんでした。

R-squaredとは、決定係数と呼ばれ、「説明変数を用いてどの程度、目的変数の分散を説明できるか」という説明力を意味しています。のR-squaredは0.8と1に近いため、まずまずのモデルと言えるでしょう。

次に、exp_lの寄与は有意ではなかったため、exp_lをはずして、exp_mを説明することを試みます。

#exp_mを目的変数、exp_n, exp_lを説明変数として重回帰分析を行う

> summary(lm(exp_m~exp_n))

Call:

lm(formula = exp_m ~ exp_n)

Residuals:

    Min 1Q Median 3Q Max 

-10.069 -8.693 -6.009 8.438 19.493 

Coefficients:

                      Estimate  Std. Error    t value Pr(>|t|)    

(Intercept)      -9.7452   17.2506     -0.565   0.587616    

exp_n              0.6113   0.1197        5.107   0.000921 ***

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 12.8 on 8 degrees of freedom

Multiple R-squared: 0.7653,Adjusted R-squared: 0.736 

F-statistic: 26.08 on 1 and 8 DF, p-value: 0.0009215 

#散布図に回帰直線を重ねて描写

> plot(exp_n, exp_m)

> abline(lm(exp_m~exp_n)

R-squredは0.06程度低下しましたが、依然として良いモデルと言えます。

実際の研究ではもっと多くの変数の組み合わせを検討してモデル探索を行い、決定係数をもとに独立変数の数を絞ります。

【参考文献】

山田剛史ほか『Rによるやさしい統計学』オーム社 2008 191 - 198pp

LaTeXItを使って簡単に数式を作成する！！！

ブログやプレゼンテーションスライドを作っているときに、一行だけ数式が欲しいということがあります。texを使うのはちょっと面倒だなあというときにおすすめなのが、LaTexItというソフトです。MacTexをインストールすると自動的にインストールされます。

MacTex
http://www.tug.org/mactex/
MacTeX.pkg 


下の画面の中で、白い部分に数式を入力して、LaTex it!をクリックすれば、数式が出力されます。「File」->「Export Image」で保存できます。

これが、LaTeXItで作成した数式です。png形式で出力しています。

対応のあるt検定とANOVA

今回は、「対応のあるt検定」について考えます。対応のあるt 検定を意識したデータ構造は以下のようになります。下の図は、ある薬剤を投与前と投与後の体重の変化を示しています。検証したい命題は「ある薬剤を投与した前後で体重は変化するか」です。

対応のないt検定と比較する目的で、下の表を提示します。各群のデータは順番は入れ替えてあるものの、上の表と同じです。薬剤投与を受けた人と薬剤投与を受けなかった人は別の人を想定します。つまり、上の表の2倍の人数が下の表では動員されていることとなります。

これから、Rを用いて対応のある場合と対応のない場合で、薬剤投与による効果

> pre <- c(67, 70, 78, 65, 75, 77, 70, 76, 79, 83)

> post <- c(66, 68, 78, 63, 74, 76, 68, 75, 78, 82)

> negative <- c(78,67,75,65,70,77,70,83,79,76)

> positive <- c(68,66,78,75,82,78,76,63,68,74)

#ドットプロット

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

#対応あり

> annot <- factor(c(rep("pre", 10), rep("post", 10), rep("negative", 10), rep("positive", 10)))

> dat <- c(pre, post, positive, negative)

> dot.plot(annot, dat)

#対応のあるt検定

> t.test(pre, post, paired=T)

 Paired t-test

data: pre and post 

t = 6, df = 9, p-value = 0.0002025

alternative hypothesis: true difference in means is not equal to 0 

95 percent confidence interval:

 0.7475686 1.6524314 

sample estimates:

mean of the differences 

                    1.2 

#対応のないt検定

> t.test(negative, positive, paired=F)

 Welch Two Sample t-test

data: negative and positive 

t = 0.4494, df = 17.91, p-value = 0.6585

alternative hypothesis: true difference in means is not equal to 0 

95 percent confidence interval:

 -4.411489 6.811489 

sample estimates:

mean of x mean of y 

     74.0 72.8 

驚くことに対応のあるt検定では有意差が検出されたのに、対応のないt検定では有意差が検出されませんでした。

この理由は、個人によるばらつきを対応のあるt検定によって制御できたことにあります。対応のないt検定の場合、個人差によるばらつきが残差に入ってしまうため、残差の分散が大きくなり、その結果薬剤の投与という要因による効果を検出するには至らなかったのです。これはまさにANOVAの考え方そのものです。というか、一般線形モデルの本質的な部分です。

対応のある一元配置分散分析を行うことで、薬剤と人という二つの要因による効果の大きさについて検証することができます。

> dat

 [1] 67 70 78 65 75 77 70 76 79 83 66 68 78 63 74 76 68 75 78 82

> drug <- factor(c(rep("negative", 10), rep("positive",10)))

> drug

 [1] negative negative negative negative negative negative negative negative

 [9] negative negative positive positive positive positive positive positive

[17] positive positive positive positive

Levels: negative positive

> person <- factor(rep(c("Hiroshi", "Takeshi", "Yamato", "Keita", "Tetsuya", "Shinta", "Yohei", "Daichi", "Yusuke", "Daisuke"),2))

> person 

[1] Hiroshi Takeshi Yamato Keita Tetsuya Shinta Yohei Daichi Yusuke 

[10] Daisuke Hiroshi Takeshi Yamato Keita Tetsuya Shinta Yohei Daichi 

[19] Yusuke Daisuke

10 Levels: Daichi Daisuke Hiroshi Keita Shinta Takeshi Tetsuya Yamato ... Yusuke

#drugとpersonの二つの要因を考慮してANOVAを行う。ただし、drugとpersonの間に交互作用はないものとする。

> summary(aov(dat~drug+person)) 

               Df  Sum Sq     Mean   Sq      F value    Pr(>F)    

drug         1           7.2              7.20           36.0    0.000202 ***

person      9        639.8            71.09         355.4    2.16e-10 ***

Residuals 9             1.8             0.20                     

---

Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

#drugだけを要因として考えてANOVA

> summary(aov(dat~drug))

               Df  Sum Sq     Mean   Sq      F value    Pr(>F)   

drug         1           7.2              7.20        0.202     0.658

Residuals  18      641.6            35.64     

drugとpersonを考慮してANOVAを行うと残差(Residuals)の平方和(Mean Sq)が小さくなり、その結果F valueが大きくなることがわかります。F valueを求めるには、

F value = (drugのMean Sq) / (ResidualsのMean Sq)

です。

このように、対応のあるt検定はanovaにおいて、要因を追加することにより系統誤差を制御していることになるのです。その結果、説明できない誤差(すなわち残差)が小さくなるため、F値は大きくなり、その結果、小さなp値が出やすくなります。

以上のことから、実験結果に影響を及ぼす要因はできるだけ見いだしておいて、制御できるものは実験デザインによって制御することが大事であると言えます。

【参考文献】
山田剛史ほか『Rによるやさしい統計学』オーム社 2008 150-156, 175-179pp

フィッシャーの三原則

統計学の父であるフィシャー：Ronald A. Fisherは、農業試験の適切な実施のための方法論として実験計画法という学問分野を確立しました。

Fisherは局所管理、反復、無作為化の３つが大切であると述べています。これはフィッシャーの３原則と呼ばれ、農学以外の幅広い分野において利用されています。
以下に、肥料の種類による作物の生育具合の評価を目的とした農場試験を例に取って説明します。

1. 局所管理（小分け）の原則
実験計画の最も原始的な形態は以下のようなものです。複数水準の間の効果を比較したいのですが、同じ場所で異なる年度で調査を行おうとしています。これでは、年度ごとの気候の差が系統誤差として入ってきます。
そこで以下のように局所管理（小分け）すると、年度ごとの気候の差を調整することができます。

2. 繰り返し（反復）の原則
小分けをしただけではまだ問題が残ります。ある水準で生育がとても良かったときに、それが果たして偶然誤差によるものなのか、水準の違い（肥料の違い）によるものなのかがはっきりしません。この問題は水準ごとに実験を反復することで、解決することができます。同じ水準内で比較することにより、偶然誤差による差を調整することができるためです。

3. 無作為化の原則
2のように水準ごとに繰り返しを作ることで、偶然誤差と系統誤差（水準によるものとそうでないものを含む）に分解することが可能となりました。最後は、実験の割り付けを無作為化することが必要となります。2.のデザインでは、水はけは日当りなど、目的要因以外の要因が、一定の方向で偏りを持って実験結果に影響を与えている可能性が高くなるからです。つまり、水準以外の系統誤差が調整できていないのです。そこで、下の図のように無作為に実験を割り付けることにより、それらの系統誤差を調整することが可能となります。

以上、フィッシャーの3原則でした。
基本的な発想は、「目的となる要因以外からの影響が極力排除すること」です。肥料以外の要因は一定になるように努力するとともに、もしそれを制御しきれないのならば、せめて肥料以外の要因が実験結果に偏って現れないようにデータ収集することが大事なのです。

この発想は臨床試験や生物学実験においても大切な発想であると思います。

【参考文献】

栗原　伸一『入門統計学』Ohmsha, 2011, 88-89, 130-138pp

二群比較と回帰分析の関係

今回は二群比較と回帰分析の関係について考えたいと思います。
前回のブログでは、ANOVAと二群比較が共に「目的とする効果による変動が誤差による変動に対して有意に大きいか」という発想を根底に持つという意味で共通しているのだという議論を行いました。
実は、回帰分析と二群比較も同様の議論が可能です。
二群比較と単回帰分析は本質的に同じことをしているのです！
単回帰分析モデルの一般式は
y = βx + α + ε
と表されます。
これは、まさに目的の効果と誤差を分離して、説明変数を説明しようとしている式です。
これを前提にして、二群比較のそれぞれの群に0と1のダミー変数を与えて、回帰分析を行うと、通常のt検定をと同じp値が得られます。
以下では、Rを用いて検証を行いました。

> cityA <- c(178, 182, 181, 179, 178, 182, 181, 179, 173, 187, 170, 190)
> cityB <- c(168, 172, 174, 166, 169, 171, 165, 175, 162, 178, 164, 176)

> dat <- c(cityA, cityB)
> dat
 [1] 178 182 181 179 178 182 181 179 173 187 170 190 168 172 174 166 169 171 165
[20] 175 162 178 164 176

#ダミー変数を用いてcityAに0を、cityBに1を与える

> city <- c(rep(0,12), rep(1,12))
> city
 [1] 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1 1 1

#この段階で散布図を描写

> plot(dat~city)

#回帰分析

> summary(lm(dat~city))

Call:

lm(formula = dat ~ city)

Residuals:

   Min     1Q Median     3Q    Max 

 -10.0   -2.5    0.0    2.5   10.0 

Coefficients:

            Estimate Std. Error t value Pr(>|t|)    

(Intercept)   180.00       1.52 118.416  < 2e-16 ***

city          -10.00       2.15  -4.652 0.000123 ***

---

Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 5.266 on 22 degrees of freedom

Multiple R-squared: 0.4959, Adjusted R-squared: 0.473 

F-statistic: 21.64 on 1 and 22 DF,  p-value: 0.0001228 

#通常の２群間比較をt検定にて行う

> t.test(cityA, cityB, var.equal=F)

Welch Two Sample t-test

data:  cityA and cityB 

t = 4.6518, df = 21.96, p-value = 0.0001233

alternative hypothesis: true difference in means is not equal to 0 

95 percent confidence interval:

  5.541324 14.458676 

sample estimates:

mean of x mean of y 

      180       170 

#散布図の上に回帰直線を描写する

> plot(city, dat)

> abline(lm(dat~city))

大学の一般教養の講義や導入書では、anovaやtテスト、回帰分析を別個の物として扱うに留まっていることが多いように感じます。しかし、一般線形モデルを意識しながらこれらを眺めると、「目的の効果と誤差を分離する」という共通原理が浮き彫りになってきます。

Peter Dalgaard『Rによる医療統計学』2007　丸善株式会社  92 - 96pp

二群比較とANOVAの関係

今回は、平均値の二群比較が、対応の無い２群の一元配置分散分析を行うことに等しいことを直感的に理解することを目指します。


今、CityAとCityBに住んでいる住人をそれぞれの群から無作為に6人だけ選出して身長を測定したとします。それが以下のデータです。
このデータに対して対応の無い２群の一元配置分散分析を行います。

分散分析とは総変動を以下のように

総変動 = 目的要因変動 + 誤差変動

と分けることが発想の基盤となります。これにしたがって、以下の図のように目的要因変動(都市による変動)と誤差による変動を計算します。

総変動を分解した結果が以下のようになります。群間変動というのは目的要因変動(都市による変動)のことであり、群変動とは誤差による変動のことを指します。

これらから、各偏差平方を計算すると以下のようになります。

分散分析とは、目的となる要因効果(群間変動)の分散が、誤差効果(標本間変動をのぞいた群内変動)の分散に比べて有意に大きいかを検定する統計手法です。これを行うにはこれら二つの比を取れば良いのです。そしてこの比は等分散の検定の時に用いた、F値です。自由度はそれぞれの自由度とします。

今回のケースで実際に計算してみると、

F = 300 / 6 = 50

となります。
自由度v1 = 1, v2 = 10の下、上側確率が5%になるF値は、,F分布表を参照すると4.96となります。
そのため、「今回のケースは都市による変動が誤差による変動を危険率5%で有意に上回る」と言うことができます。

一方、今回のケースにStuden t検定を行うとどうなるでしょうか。個々から下はRを用いて計算を行います。

$ R
> cityA <- c(178, 182, 181, 179, 178, 182)
> cityB <- c(168, 172, 174, 166, 169, 171)

> t.test(cityA, cityB, var.equal=T)

Two Sample t-test

data:  cityA and cityB 

t = 7.0711, df = 10, p-value = 3.411e-05

alternative hypothesis: true difference in means is not equal to 0 

95 percent confidence interval:

  6.848936 13.151064 

sample estimates:

mean of x mean of y 

      180       170 

p< 0.05なので、二群の平均値の間には有意な差があると言うことができます。

ところで、t値は7.0711と求まりました。

実はこのt値を2乗してみると面白い結果を得ることができます。

> 7.0711 * 7.0711

[1] 50.00046

なんと驚いたことに、先ほど計算したF値とStudent t testで計算したt値の二乗が一致するではありませんか。これはただの偶然ではないのです。

一般に、平均値の二群比較でもとめられるt値の二乗は、対応の無い２群の一元配置分散分析を行うF値に等しいのです。

このことから、tテストを行うということはANOVAを特定下で行うことに等しいということができます。

この話題は、結局平均値の二群比較とANOVAと回帰分析を巻き込んで、一般化線形モデルというトッピックに発展するようです。

まだ勉強の途中でよく理解できていないので、あまり細かいことは良く理解できていません。しかし、現状理解できているのは、「見たい要因による効果」と「偶然による誤差」を一次結合の形式で分離することができるというモデル(信念)を頭において、これらの理論は構築されているということです。
非線形であった場合、例えば誤差がある値を超えてきた場合に、見たい要因の変動にもとても大きな(小さな)な影響を与えるようなことが起きると思われます(一例ですが)。非線形の場合は、さまざまな要因間が絡み合っていると考えるはずです。要するに、「世の中、物事の足し算でできているような単純なものではない」とすることではないでしょうか。

少し、哲学的になってきてしまいましたが、現状私が理解できている(しているつもり)の範囲を書きました。

【参考文献】

栗原　伸一『入門統計学』Ohmsha, 2011, 88-89, 130-138pp

霊長類の進化

先週、愛知県犬山市にある日本モンキーセンター(http://www.j-monkey.jp/)を見学してきました。様々な種類の霊長目の動物達が活発に動き回る姿を目の当たりにして、私たちヒトという動物の生きる意味について問い直すきっかけとなりました。

さて、今回は霊長目(霊長類)の進化の歴史について勉強したことをまとめます。

例調理は、中生代の白亜紀末期、今からおよそ7000万年前にほかのグループから分かれたと考えられています。
霊長類の独自の特徴は以下のようなものが挙げられます。

• 樹上への適応を通じた独自の進化
• 平爪を有する
• 発達した鎖骨
• 母指対立により物を握ることが可能
• 両眼で前方の同じものを見ることが可能
• 赤色オプシン遺伝子を重複させて変異を加えることで、緑色オプシンを復活させている

原始的な真アルコンタ類からヒトに至るまでの進化の過程で重要な分岐は以下の通りです。

• 約7000万年前(白亜紀終わり頃)
• 霊長目がツパイやヒヨケザルの仲間から分岐
• 約3000万年前(新生代)
• 類人猿が旧世界ザルと分岐
• 約700万年前
• ヒトとチンパンジーが分かれた
• 約300万年前
• チンパンジーとボノボが分岐

日本モンキーセンターのゴリラやチンパンジーの檻のそばにはヒト様の檻も用意されていて、中にはベンチが置いてあり記念撮影ができるようになてちます。その檻に張られたヒトの紹介文の内容は、ヒトの可能性や危険性に言及したもので大変示唆に富んでいました。よって、ここではこれを引用させて頂きます。

ヒト

英名　：　Human

学名　：　Homo sapiens

分布　：　地球上の陸地のほぼ全域

身体

• しっぽがなく、からだの毛は少ない。
• 器用な前足と大きな脳を持つ。
• 二本足でまっすぐ立って歩き、サルの仲間のくせに木登りがへた。

食性

• なんでもよく食べる雑食性だが、１頭ごとの好き嫌いは多い。
• 食べ物をいろいろ加工したり、貯めこんだりする。

習性

• 好奇心が強く、遊び好き。危険な眼にあっても、なかなかこりない。
• 大きく複雑な群れを作って生活し、音声によるコミュニケーションがさかんである。
• 「文字」や「絵」などを通じて情報をやりとりし、いつも仲間とふれ合うことをもとめている。
• しばしば、いさかいをおこあし、仲間同士で殺し合いをすることも珍しくない。
• 本能だけではうまく生きられず大人になるまでに身につけなくてはならないことがとげも多い。
• 優れた知能をもち、大きな可能性をひめているるが、失敗したときにまわりに与える影響も大きい。それを忘れたとき地球上でもっとも危険な動物となる。

【参考文献】

井出利憲「分子生物学講義中継 生物の多様性と進化の驚異」羊土社 2010 310-311pp

財団法人日本モンキーセンター「霊長類　ガイドブック」 2012 64pp

Rを用いた単回帰分析の基本

ある試薬を培養細胞に添加したときのある遺伝子の発現量を制御するための実験を行った。実験結果は以下の通りである。このデータを用いて単回帰分析を行う。

以下解析。

#Rを起動
$ R

#データの入力
x <- c(0.1, 0.2,0.3,0.4,0.5,0.6)
y<-c(3.1,3.8,5.7,6.7,7.5,9.3)
dat <- data.frame(x, y)

#データの確認
dat

    x   y
1 0.1 3.1
2 0.2 3.8
3 0.3 5.7
4 0.4 6.7
5 0.5 7.5
6 0.6 9.3

#データの描写

png("dat1.png"); plot(dat, xlab="reagent dose(ug)", ylab="gene expresion"); dev.off();

#回帰分析

model <- lm(y~x, data=dat)

#回帰分析結果の要約表示

summary(model)

Call:

lm(formula = y ~ x, data = dat)

Residuals:  #各残差

            1             2             3            4              5             6 

 0.16190 -0.36952  0.29905  0.06762 -0.36381  0.20476 

Coefficients:

                 Estimate     Std. Error    t value    Pr(>|t|)            #推定値

(Intercept)   1.7067        0.3056      5.585     0.00504 **     #定数項(β0)

x                12.3143       0.7847    15.693     9.63e-05 ***  #回帰係数(β1)

---

Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 0.3283 on 4 degrees of freedom  #残差の標準誤差

Multiple R-squared: 0.984, Adjusted R-squared:  0.98 #寄与率、調整済寄与率

F-statistic: 246.3 on 1 and 4 DF,  p-value: 9.632e-05        #F0値、p値

#回帰係数の95%信頼区間を算出

confint(model, level=0.95)

                          2.5 %      97.5 %
(Intercept)  0.8581746    2.555159    #β0の95%信頼区間の下限と上限
x              10.1355593   14.493012   #β1の95%信頼区間の下限と上限

#実験データの散布図と回帰直線の描写

png("dat2.png");

plot(dat, xlab="reagent dose(ug)", ylab="gene expresion"); 

abline(model, col="red")

dev.off();

【参考文献】

辻谷 将明、 和田 武夫『Rで学ぶ確率・統計』共立出版 2012 122 -136 pp.