データ処理初めの1歩(50分) - ryamadaの遺伝学・遺伝統計学メモ

大学１回生に50分でデータ処理の「実務」を教えるとき
コンテンツ
- 記録(スプレッドシート)
  - {0,1}
  - 連続値
- Rへの読み込み・取り込み
- 集計・記述統計
- 視覚化
- 検定
記録(スプレッドシート、エクセル)
- ファイル(アクセスしてダウンロード)
- テキスト形式で保存しよう

- - テキスト形式ファイル(アクセスしてダウンロードしてもよい)
- 保存したら中身を確認しよう

Rへの読み込み → 読み込みができるようになるために(脇道。必要に応じて読む)
そもそもRの起動と終了
- 今日は、それをスキップして、読み込んでしまう

#data1 <- read.table("data1.txt",header=TRUE)
data1 <- read.table("http://www.statgenet.med.kyoto-u.ac.jp/StatGenet/lectures/2016/MdPhDStatPrimer/Data1.txt",header=TRUE)

中身を見てみよう

data1

分割表にしてみる

data1.table <- table(tr=data1$Treatment,re=data1$Reaction)
data1.table

分割表は検定してみる

chisq.test(data1.table)
fisher.test(data1.table)

連続値データ
- ファイル(アクセスしてダウンロード)
  - テキスト形式ファイル
- エクセルファイルをダウンロードして自分でテキストファイル化する練習をしてもよい

data2 <- read.table("http://www.statgenet.med.kyoto-u.ac.jp/StatGenet/lectures/2016/MdPhDStatPrimer/Data2.txt",header=TRUE)

中身を見てみよう

data2

箱ひげ図を描いてみる

boxplot(data2$Reaction ~ data2$Treatment)

２群を別々のオブジェクトに格納してみる

data2.0 <- data2$Reaction[which(data2$Treatment==0)]
data2.1 <- data2$Reaction[which(data2$Treatment==1)]

群ごとに値の様子を見る

hist(data2.0)
mean(data2.0)
min(data2.0)
max(data2.0)
median(data2.0)

summary(data2.0)

検定してみる

t.test(data2.0,data2.1)

データを作って、回帰してみる

# サンプル数 n
n <- 100
# 遺伝子Xの発現量は、平均120、分散15だと言う
# 乱数でデータ生成する
m <- 120
v <- 15
x <- rnorm(n, m, sqrt(v))
hist(x)

# 遺伝子Yの発現量は、y = 2.4 x + 30の関係がありつつ、正規乱数誤差(平均0、分散 5)があると言う
a <- 2.4
b <- 30

m.r <- 0
v.r <- 5
y <- a * x + b + rnorm(n, m.r, v.r)
plot(x,y)

回帰係数を計算してみる

cor(x,y)

線形回帰をしてみる

lm.out <- lm(y ~ x)
lm.out
lm.out$coefficients
lm.coef <- lm.out$coefficients
plot(x,y)
abline(lm.coef[1],lm.coef[2],col="red")
points(x,predict(lm.out),pch=20,col="green")

この後、どうする