Hazy Ideas

日々の勉強の気づきを書き出しています

【Rで疫学】記述統計表とグラフ

疫学 R

今回は、適当に作成したデータを用いて、記述統計に用いる要約統計量をまとめた表と、可視化の方法として箱ひげ図を作成してみます。

方針としては、子供の肺機能と身長・体重、既往歴などを含めたデータを作成して、学年と性別ごとに並べて見比べることにします。使用するＲパッケージは、"モダンな"データ形成ができるtidyverseと、直感的な作図ができるesquisseです。

現実ではきれいなデータはあり得ませんが、今回はテストなので性比や身長、FVCの値は適当に指定していきます。

データ形成

#libraryインポート
library(tidyverse)
library(esquisse)

#乱数種の指定
set.seed(1)

#条件
#中学１～３年生、100人ずつ
#性別 male/female


#身長 130～180cmくらい
#体重 40～60kgくらい
#肺機能として努力性肺活量（FVC）3.5～5.5Lくらい
#病気Aの既往歴　yes/no

#1年生 100人
id1 <- c(1:100)
grad1 <- rep(1,100)
#今回は性別はランダムに割り振ります（男女比半々）
s <- c("male", "female")
sex1 <- sample(s, 100, replace=T)
#平均身長120cm、標準偏差20とします
height1 <- rnorm(n=100, mean =120, sd=20)
#平均体重45kg、標準偏差10とします
weight1 <- rnorm(n=100, mean =45, sd=10)
#平均2.5L、標準偏差0.3とします
FVC1 <- rnorm(n=100, mean =2.5, sd=0.3)
#病気Aの既往歴ありの割合は20％とします
a <- c("yes", "no")
sickA1 <- sample(a, 100, replace=T, prob=c(0.2, 0.8)) 

G1 <- data.frame(id1, grad1, sex1, height1, weight1, FVC1, sickA1)

sample()関数は、乱数を発生させランダムな割付を行う関数です。結果が毎回変わるため、今回はset.seed(1)を指定して、結果の再現性があるようにしました。他2学年については以下のようにしました。

#2年生 100人
id2 <- c(101:200)
grad2 <- rep(2,100)
#今回は性別はランダムに割り振ります（男女比半々）
sex2 <- sample(s, 100, replace=T)
#平均身長135cm、標準偏差20とします
height2 <- rnorm(n=100, mean =135, sd=20)
#平均体重48kg、標準偏差10とします
weight2 <- rnorm(n=100, mean =48, sd=10)
#平均3.3L、標準偏差0.4とします
FVC2 <- rnorm(n=100, mean =3.3, sd=0.4)
#病気Aの既往歴ありの割合は25％とします
sickA2 <- sample(a, 100, replace=T, prob=c(0.25, 0.8)) 

G2 <- data.frame(id2, grad2, sex2, height2, weight2, FVC2, sickA2)

#3年生 100人
id3 <- c(201:300)
grad3 <- rep(3,100)
#今回は性別はランダムに割り振ります（男女比半々）
sex3 <- sample(s, 100, replace=T)
#平均身長150cm、標準偏差20とします
height3 <- rnorm(n=100, mean =150, sd=20)
#平均体重52kg、標準偏差10とします
weight3 <- rnorm(n=100, mean =52, sd=10)
#平均4.0L、標準偏差0.5とします
FVC3 <- rnorm(n=100, mean =4.0, sd=0.5)
#病気Aの既往歴ありの割合は30％とします
sickA3 <- sample(a, 100, replace=T, prob=c(0.30, 0.8)) 

G3 <- data.frame(id3, grad3, sex3, height3, weight3, FVC3, sickA3)

最後にこれらのデータをまとめて一つのデータフレームにします。その際に、列名が異なると結合できないため、同名に変更しました。

colnames(G1) <- c("ID","grade","sex","height","weight","FVC","medhis")
colnames(G2) <- c("ID","grade","sex","height","weight","FVC","medhis")
colnames(G3) <- c("ID","grade","sex","height","weight","FVC","medhis")

DF <- rbind(G1, G2, G3)

これで仮データの作成が完了しました。いよいよ次に、記述統計表（要約統計量の表）を作成します。

記述統計表の作成

まずはデータの全体像を眺めてみます。見やすくはないですが、基本のベースRのplot()を使います。

#全体像の確認
plot(DF)

学年（2列目）を見ると、学年が上がるごとに身長・体重・FVCが上がっていることが読み取れます（そういう風にデータ形成したので当たり前ですが）。

次にデータの数値的な様子を読み解くために、統計量を出してみます。まずはベースRのsummary()を使います。

#全体の統計値
summary(DF)
       ID              学年       性別          
 Min.   :  1.00   Min.   :1   Length:300        
 1st Qu.: 75.75   1st Qu.:1   Class :character  
 Median :150.50   Median :2   Mode  :character  
 Mean   :150.50   Mean   :2                     
 3rd Qu.:225.25   3rd Qu.:3                     
 Max.   :300.00   Max.   :3                     
      身長             体重            FVC       
 Min.   : 81.71   Min.   :15.87   Min.   :1.800  
 1st Qu.:115.70   1st Qu.:41.68   1st Qu.:2.700  
 Median :134.16   Median :48.62   Median :3.243  
 Mean   :134.88   Mean   :48.26   Mean   :3.275  
 3rd Qu.:151.95   3rd Qu.:55.70   3rd Qu.:3.766  
 Max.   :193.79   Max.   :74.06   Max.   :5.115  
    既往歴         
 Length:300        
 Class :character  
 Mode  :character

これだけ見ても仕方ないですよね・・・。学年や性別などで分けてみたいと思います。パッケージtidyverseに含まれている、ggplot()を用います。パイプ演算子%>%や、使われている関数（dplyrの関数）の詳細説明は、他に譲らせていただきます。

summary <- 
  DF %>%
  group_by(grade) %>%
  summarize(
    No. = n(), #総数
    Male = sum(sex=="male"), #男女比として男性を数える
    Height = mean(height),　#平均身長
    Weight =mean(weight),　#平均体重
    FVC = mean(FVC),　#平均FVC
    MedHis = sum(medhis=="yes")　#病歴ありの数
  )
summary

学年ごとの性別比（男性数）、身長、体重、FVCの平均値が一つの表に表示されました。

グラフの作成

では次にグラフ化してみます。今回、学年が数字numericで、このままではグラフにできないため文字列characterへ事前に変換しておきます。

#numericからcharacterへ変換
DF2 <- DF %>% mutate(grade = as.character(grade))

グラフ作成のポイントですが、どんな図を作りたいかあらかじめ考えておくことです。それを直感的操作で手助けしてくれるのが、esquisse()です。動かしてみる便利さがわかります。

#作図できるパターンを直感的に理解する
esquisser(DF2)

上記の画面が表示されるので、ｘ軸にgrade、ｙ軸にweightをドラッグ＆ドロップします。右下のCode▲を押すと、表示されているグラフのコードが表示されます。

このコードをそのまま使ってもいいですし、これをベースに書き換えていくのでもよいです。ではまずは学年別身長をグラフにしてみます。

#学年別身長
g_height <- ggplot(DF2)+
  geom_boxplot(aes(x=grade, y=height))
g_height

次に学年別体重を、男女別に表示してみます。

#学年別体重、男女別
g_weight <- ggplot(DF2)+
  geom_boxplot(aes(x=grade, y=weight, fill=sex))
g_weight

では最後に、学年別FVCを、男女別、かつ既往歴有無をグループ分けして表示します。

#学年別FVC,男女別,既往歴で分ける
g_FVC <- ggplot(DF2)+
  geom_boxplot(aes(x=grade, y=FVC, fill=sex))+
  facet_wrap(vars(medhis))
g_FVC

まとめ

今回は記述統計表とグラフ作成の練習を紹介しました。今回はデータを適当に作成したため、性別や既往歴で違いは見れませんでしたが、図表作成の練習にはなるかと思いますので、ぜひ練習してみてください。