J

ggglot2 のさまざまなグラフの例（統計的観点から）

R ggplot2

f:id:joker1110:20120215053727p:image:w360

背景と目的

さまざまな統計データを扱う上で、グラフなどへの可視化はとても重要な手段である。ggplot2 は統計的ツールとして有名な R のグラフライブラリであり、本ブログでもその使い方を紹介している。ここでは、その具体的な例として、三重大の奥村先生が公開しているグラフの描き方に描かれている素晴らしい統計的グラフの数々を、練習をかねて ggplot2 で描いてみる。

テーマのカスタマイズ

ggplot2 ではデフォルトでも十分美しいグラフが描ける。ここでは、グラフの描き方により近いグラフを描くために、一貫してggplot2 のテーマを作成 - joker8phoenix の日記のテーマを利用する。このテーマをmyggplot2.Rとして保存して、

source("~/[ファイルの場所]/myggplot2.R")

とすれば、mytheme()というテーマを定義できる。

棒グラフ・柱状グラフ

地域別面積
  • データ
region = c("北海道","本州","四国","九州","沖縄")
area = c(83457,231113,18792,42191,2276)/10000
  • グラフ

f:id:joker1110:20120212030156p:image:w360

  • ソース
source("~/[ファイルの場所]/myggplot2.R") # 自分の環境に合わせる
quartz(width=8.267717,height=5.826772) # Macのみ利用可能
areadata = data.frame(region=region,area=area,ord=c(1:5))
garea <- ggplot(areadata)
garea <- garea + geom_bar(aes(x=reorder(region,ord),weight=area),fill="grey",colour="black")
garea <- garea + labs(x=NULL, y=expression(paste("面積（万", km^2, "）")))
garea <- garea + mytheme()
garea <- garea + coord_cartesian(ylim=c(0,25))
garea
grid.gedit("text", gp=gpar(fontfamily="HiraMaruProN-W4"))
都道府県別人口
  • データ
kenmei = c("北海道", "青森県", "岩手県", "宮城県", "秋田県", "山形県", "福島県", "茨城県", "栃木県", "群馬県", "埼玉県", "千葉県", "東京都", "神奈川県", "新潟県", "富山県", "石川県", "福井県", "山梨県", "長野県", "岐阜県", "静岡県", "愛知県", "三重県", "滋賀県", "京都府", "大阪府", "兵庫県", "奈良県", "和歌山県", "鳥取県", "島根県", "岡山県", "広島県", "山口県", "徳島県", "香川県", "愛媛県", "高知県", "福岡県", "佐賀県", "長崎県", "熊本県", "大分県", "宮崎県", "鹿児島県", "沖縄県")
population = c(5535, 1392, 1352, 2340, 1108, 1188, 2052, 2964, 2011, 2012, 7113, 6122, 12838, 8917, 2391, 1101, 1168, 812, 871, 2171, 2100, 3800, 7403, 1875, 1402, 2629, 8806, 5586, 1404, 1012, 595, 725, 1948, 2869, 1463, 794, 1003, 1444, 773, 5054, 856, 1440, 1821, 1200, 1136, 1717, 1376)
  • グラフ

f:id:joker1110:20120212031122p:image:w360

  • ソース
quartz(width=8.267717,height=5.826772)
pdata = data.frame(kenmei=kenmei,population=population)
gp <- ggplot(pdata)
gp <- gp + geom_bar(aes(x=reorder(kenmei,population),weight=population/1000),fill="grey",colour="black")
gp <- gp + scale_y_continuous(breaks=population)
gp <- gp + labs(y=NULL, x=NULL)
p <- seq(0,12,2)
s = sum(population) / 1000
t = seq(0,10,2)
gp <- gp + geom_linerange(aes(x=0,ymin=0,ymax=max(population)/1000),colour="grey50")
gp <- gp + geom_segment(aes(x=0,y=p,xend=-1,yend=p))
gp <- gp + geom_text(aes(x=0,y=p,label=p),vjust=2)
gp <- gp + geom_linerange(aes(x=48,ymin=0,ymax=s*max(t)/100),colour="grey50")
gp <- gp + geom_segment(aes(x=48,y=s*t/100,xend=49,yend=s*t/100))
gp <- gp + geom_text(aes(x=48,y=s*t/100,label=paste(t,"%",sep="")),vjust=-1,hjust=0.2)
gp <- gp + geom_text(aes(x=-6,y=s*max(t)/100/2,label="人口（百万人）"))
gp <- gp + geom_text(aes(x=54,y=s*max(t)/100/2,label="人口の割合"))
gp <- gp + mytheme(base_size=8)
gp <- gp + opts(axis.line=theme_blank())
gp <- gp + coord_flip(ylim=c(-0.2,14))
print(gp)
grid.gedit("text", gp=gpar(fontfamily="HiraMaruProN-W4"))
国立大学運営費交付金
  • データ
year = 2004:2011
money = c(12415, 12317, 12214, 12043, 11813, 11695, 11585, 11528)
  • グラフ

f:id:joker1110:20120212032157p:image:w360

  • ソース
quartz(width=8.267717,height=5.826772)
mdata = data.frame(year=year,money=money)
gm <- ggplot(mdata)
gm <- gm + geom_bar(aes(x=as.character(year),y=money),fill="grey",colour="black")
gm <- gm + geom_text(aes(x=as.character(year),y=money,label=money), vjust=-0.5)
gm <- gm + mytheme(title="国立大学法人運営費交付金予算",base_size=14)
gm <- gm + labs(x=NULL,y="億円")
gm <- gm + coord_cartesian(ylim=c(0,14000))
gm
grid.gedit("text", gp=gpar(fontfamily="HiraMaruProN-W4"))
  • 折れ線バージョン

f:id:joker1110:20120212032651p:image:w360

  • ソース
quartz(width=8.267717,height=5.826772)
gm <- ggplot(mdata)
gm <- gm + geom_point(aes(x=year,y=money),size=3)
gm <- gm + geom_line(aes(x=year,y=money))
gm <- gm + mytheme(title="国立大学法人運営費交付金予算")
gm <- gm + labs(x=NULL,y="億円")
gm
grid.gedit("text", gp=gpar(fontfamily="HiraMaruProN-W4"))
交通事故死者数
  • データ
year = 1995:2009
death = c(10679, 9942, 9640, 9211, 9006, 9066, 8747, 8326, 7702, 7358, 6871, 6352, 5744, 5155, 4914)
  • グラフ

f:id:joker1110:20120214204658p:image:w360

  • ソース
quartz(width=8.267717,height=5.826772)
ddata <- data.frame(year=year,death=death)
gd <- ggplot(ddata)
gd <- gd + geom_bar(aes(x=as.character(year),y=death),fill=cCyan,colour="black")
gd <- gd + scale_x_discrete(breaks=year,labels=substr(year,3,4))
gd <- gd + labs(x=NULL,y="人")
gd <- gd + mytheme(base_size=20,margin.top=3,margin.side=1,title="交通事故死者数")
gd <- gd + coord_cartesian(ylim=c(0,12000))
gd
  • 折れ線バージョン

f:id:joker1110:20120214204655p:image:w360

  • ソース
quartz(width=8.267717,height=5.826772)
gd <- ggplot(ddata)
gd <- gd + geom_line(data=data.frame(year=c(1970,1995),death=c(16765,10679)),aes(x=year,y=death),colour="grey",linetype=2)
gd <- gd + geom_point(data=data.frame(year=1970,death=16765),aes(x=year,y=death),size=4)
gd <- gd + geom_point(aes(x=year,y=death),size=4)
gd <- gd + geom_line(aes(x=year,y=death))
gd <- gd + mytheme(20)
gd <- gd + opts(title="交通事故死者数")
gd <- gd + labs(x=NULL,y="人")
gd <- gd + scale_x_continuous(breaks=c(1970, 1995, 2000, 2005, 2009))
print(gd)
所得の分布
  • データ
income=c(6.2,11.7,12.9,13.2,10.6,9.2,7.8,6.3,4.9,4.2,2.9,2.2,1.9,1.3,1,0.7,0.6,0.4,0.4,0.2,1.4)
  • グラフ

f:id:joker1110:20120215043551p:image:w360

  • ソース
t=c(0:20)
idata <- data.frame(income=income,t=t,it=as.character(t))
quartz(width=8.267717,height=5.826772)
gi <- ggplot(idata)
gi <- gi + geom_rect(aes(xmin=t,xmax=t+1,ymin=0,ymax=income),colour="black",fill="grey90")
gi <- gi + geom_text(aes(x=t,y=income, label=paste(income,"%",sep="")),vjust=-0.4,hjust=-0.2,size=4)
gi <- gi + geom_text(aes(x=t, y=0, label=paste(seq(0,2000,100),"万",sep="")),hjust=1.1,size=4)
gi <- gi + geom_linerange(aes(x=5.668,ymin=0,ymax=1))
gi <- gi + geom_text(aes(x=5.668, y=1, label="平均値566.8万円"),hjust=0,size=3)
gi <- gi + geom_linerange(aes(x=4.51,ymin=0,ymax=1))
gi <- gi + geom_text(aes(x=4.51, y=1, label="中央値451万円"),hjust=0,size=3)
gi <- gi + scale_x_discrete(breaks=t)
gi <- gi + scale_y_continuous(breaks=NULL)
gi <- gi + mytheme()
gi <- gi + labs(y=NULL, x=NULL)
gi <- gi + opts(title="2006年 所得分布",axis.line=theme_blank(),axis.ticks=theme_blank(),axis.text.x=theme_blank())
gi <- gi + coord_flip(ylim=c(-2,16))
gi
grid.gedit("text", gp=gpar(fontfamily="HiraKaku"))
賃金の分布
  • データ
salary = c(0,100,120,140,160,180,200,220,240,260,280,300,
         320,340,360,400,450,500,600,700,800,900,1000,1200)
rate = c(0.2,0.5,1.3,2.6,4.4,5.6,7.5,7.6,7.6,6.9,6.4,6.2,
         5.4,4.8,8.1,7.6,5.2,6.2,3.0,1.4,0.7,0.3,0.3,0.2)
range = c(salary[2:24],1400)-salary
  • グラフ

f:id:joker1110:20120217031310p:image:w360

  • ソース
quartz(width=8.267717,height=5.826772)
sdata <- data.frame(salary=salary, rate=rate, range=range)
t = c(0,100,200,300,360,400,450,500,600,700,800,900,1000,1200)
gs <- ggplot(sdata)
gs <- gs + geom_rect(aes(xmin=salary,xmax=salary+range,ymin=0,ymax=rate/range),fill="grey",colour="black")
gs <- gs + geom_text(aes(x=salary+range/2,y=rate/range,label=paste(rate,"%","")),vjust=0.4, hjust=-0.3,size=2)
gs <- gs + scale_x_continuous(breaks=t,labels=paste(t/10,"万",sep=""))
gs <- gs + mytheme()
gs <- gs + opts(axis.line=theme_blank(),axis.ticks=theme_blank(),axis.text.x=theme_blank())
gs <- gs + labs(y=NULL, x=NULL)
gs <- gs + coord_flip(ylim=c(0,0.5))
print(gs)
grid.gedit("text", gp=gpar(fontfamily="HiraKaku"))
アンケート
  • データ
select = c("千倍","一万倍","十万倍","百万倍","千万倍","わからない")
value = c(18,12,21,21,17,12)
  • グラフ

f:id:joker1110:20120217031751p:image:w360

  • ソース
t <- c(6:1)
tdata <- data.frame(select=select, value=value, t=t)
quartz(width=8.267717,height=5.826772)
gt <- ggplot(tdata)
gt <- gt + geom_bar(aes(x=as.character(t),y=value),colour="black",fill="grey")
gt <- gt + geom_text(aes(x=t,y=value,label=paste(value,"%",sep="")),hjust=-0.5)
gt <- gt + mytheme()
gt <- gt + opts(title="兆は百万の何倍？",axis.line=theme_blank(),axis.ticks=theme_blank(),axis.text.x=theme_blank())
gt <- gt + scale_x_discrete(breaks=t,labels=select)
gt <- gt + labs(x=NULL, y=NULL)
gt <- gt + coord_flip(ylim=c(0,25))
print(gt)
アンケート（複数回答）
  • データ
製品サービス=c("パソコン","ファッション","携帯音楽プレーヤー","通信機器","国内旅行","音楽","書籍","アニメ・漫画","ゲーム","外食・食べ歩き","映画","カメラ","海外旅行","テレビ","語学、資格試験","化粧品、エステ","自動車")
割合=c(62.1, 53.9, 50.6, 49.9, 44.0, 43.7, 42.9, 42.0, 38.4, 37.6, 35.5, 35.0, 33.9, 28.8, 27.3, 26.2, 22.8)
  • グラフ

f:id:joker1110:20120217032300p:image:w360

  • ソース
pdata <- data.frame(product=製品サービス,rate=割合)
quartz(width=8.267717,height=5.826772)
gp <- ggplot(pdata)
gp <- gp + geom_bar(aes(x=reorder(product,rate),y=rate),fill="grey",colour="black")
gp <- gp + geom_text(aes(y=rate, x=product, label=paste(割合,"%",sep="")),hjust=-0.2)
gp <- gp + mytheme(title="大学生の興味関心がある製品・サービスランキング")
gp <- gp + opts(axis.line=theme_blank(),axis.text.x=theme_blank(),axis.ticks=theme_blank())
gp <- gp + labs(x=NULL,y=NULL)
gp <- gp + coord_flip(ylim=c(0,75))
gp
積み重ね棒グラフ
  • データ

ここで扱うデータはグラフの例：積み重ね棒グラフのデータをcsv形式でDataフォルダに保存してあるとする。

X = read.csv(file("Data/highschool_info.csv",encoding="cp932"),header=T)
X = cbind(X, data.frame(教えていない=100-X[,2]-X[,3], ord=X[,2]+X[,3]))
  • グラフ

f:id:joker1110:20120217033049p:image:w360

  • ソース
X.melt <- melt.data.frame(X,id.vars=c("X","ord"),variable_name="category",na.rm=T)
quartz(width=8.267717,height=5.826772)
gx <- ggplot(X.melt)
gx <- gx + geom_bar(aes(x=reorder(X,ord),y=value,fill=category),colour="black",stat="identity")
z = o[length(o)]
gx <- gx + geom_text(aes(y=X.melt[z,4]/2, label="教えている", x=25))
gx <- gx + geom_text(aes(label="一部教えている", y=X.melt[z,4]+X.melt[z+length(o),4]/2, x=25))
gx <- gx + labs(x=NULL,y=NULL)
gx <- gx + mytheme()
gx <- gx + opts(legend.position="none",axis.line=theme_blank(),axis.text.x=theme_blank(),axis.ticks=theme_blank())
gx <- gx + scale_fill_grey()
gx <- gx + coord_flip(ylim=c(0,100))
gx
grid.gedit("text", gp=gpar(fontfamily="HiraKaku"))
積み重ね棒グラフ2
  • データ
A = c(55.5, 17, 9.3, 5.7)
B = c(6.4, 4.9, 7.4, 2.6)
C = c(9.2, 4.8, 10, 3.5)
  • グラフ

f:id:joker1110:20120217033818p:image:w360

  • ソース
quartz(width=8.267717,height=5.826772)
edata <- data.frame(id=rep(1:4,1),A=A,B=B,C=C)
edata.melt <- melt.data.frame(edata,id.vars="id",variable_name="category")
ge <- ggplot(edata.melt)
ge <- ge + geom_bar(aes(x=category,weight=value,fill=factor(id)),colour="black")
gval <- cumsum(A)-A/2
ge <- ge + geom_text(aes(x="A",y=gval,label=c("1年次", "2年次", "3年次", "単位制"),colour=factor(id)))
ge <- ge + scale_colour_grey(start=0.8,end=0.2) 
ge <- ge + scale_fill_grey()
ge <- ge + scale_x_discrete(breaks=c("A","B","C"),labels=c("情報A","情報B","情報C"))
ge <- ge + labs(x=NULL,y=NULL)
ge <- ge + mytheme()
ge <- ge + opts(legend.position="none")
ge <- ge + coord_cartesian(ylim=c(0,100))
ge
grid.gedit("text", gp=gpar(fontfamily="HiraKaku"))
迷惑メールの数
  • データ

ここで扱うデータはグラフの例：迷惑メールの数のデータをDataフォルダに保存してあるとする。

chainletters = read.csv("Data/chainletters.csv", as.is=TRUE,
                        fileEncoding="SJIS", comment.char="#")
t = as.POSIXct(chainletters$期間)
  • グラフ

f:id:joker1110:20120217034316p:image:w360

  • ソース
quartz(width=8.267717,height=5.826772)
chainletters$期間 = as.POSIXct(chainletters$期間)
gc <- ggplot(chainletters)
gc <- gc + geom_bar(aes(x=期間+3.5*24*60*60,weight=件数),binwidth=7*24*60*60,fill=cCyan,colour="white")
gc <- gc + labs(x=NULL,y=NULL)
gc <- gc + mytheme()
r <- range(chainletters$期間)
gc <- gc + scale_x_datetime(major="2 months", format="%m月",limits=c(r[1],r[2]+3.5*24*60*60))
gc <- gc + opts(title="...チェーンメールの転送数")
gc <- gc + coord_cartesian(ylim=c(0,12000))
gc

折れ線グラフ

地球温暖化
  • データ
気温 = c(-0.49, -0.56, -0.6, -0.56, -0.53, -0.32, -0.34, -0.51, -0.41, -0.34, -0.43, -0.55, -0.62, -0.69, -0.55, -0.45, -0.63, -0.67, -0.67, -0.64, -0.67, -0.58, -0.55, -0.38, -0.28, -0.49, -0.56, -0.4, -0.43, -0.36, -0.28, -0.41, -0.39, -0.41, -0.31, -0.22, -0.32, -0.33, -0.45, -0.23, -0.2, -0.23, -0.38, -0.23, -0.31, -0.33, -0.22, -0.19, -0.22, -0.17, -0.11, -0.12, -0.08, 0.03, -0.11, -0.26, -0.28, -0.26, -0.27, -0.34, -0.21, -0.15, -0.08, -0.31, -0.32, -0.41, -0.14, -0.08, -0.14, -0.18, -0.09, -0.07, -0.04, -0.34, -0.28, -0.21, -0.21, -0.23, -0.12, -0.15, -0.26, -0.14, -0.02, -0.29, -0.24, -0.33, -0.04, -0.13, -0.01, 0.02, 0.06, -0.06, 0.09, -0.09, -0.11, -0.02, 0.14, 0.12, 0.05, 0.19, 0.12, -0.02, 0.01, 0.08, 0.16, 0.06, 0.24, 0.37, 0.15, 0.15, 0.27, 0.31, 0.31, 0.27, 0.32, 0.31, 0.28, 0.2, 0.31)
  • グラフ

f:id:joker1110:20120221140844p:image:w360

  • ソース
quartz(width=11, height=3.3)
tdata <- data.frame(time=1891:2009,Temp=気温)
gt <- ggplot(tdata)
gt <- gt + geom_line(aes(x=time,y=Temp))
gt <- gt + geom_point(aes(x=time,y=Temp))
gt <- gt + mytheme(title="地球の平均気温の変化（単位：℃，相対値）")
gt <- gt + opts(panel.border=theme_rect(colour="black"))
gt <- gt + labs(x=NULL,y=NULL)
gt <- gt + scale_x_continuous(breaks=c(1891, seq(1900,2000,20), 2009))
gt
男女別睡眠時間の推移
  • データ
睡眠時間 = data.frame(
  年=seq(1986,2006,5),
  男=c(7+56/60,7+50/60,7+52/60,7+49/60,7+47/60),
  女=c(7+39/60,7+34/60,7+36/60,7+35/60,7+32/60))
  • グラフ

f:id:joker1110:20120221141447p:image:w360

  • ソース
quartz(width=4.245817,height=5.926181)                      # Mac
sleep <- melt.data.frame(睡眠時間, id.vars="年",variable_name="sex")
gs <- ggplot(sleep,aes(x=年,y=value))
gs <- gs + geom_line(aes(colour=sex)) + geom_point(aes(colour=sex))
gs <- gs + scale_colour_manual(aes(sex), value=qcolours)
gs <- gs + mytheme(title="15歳以上の日本人の睡眠時間")
gs <- gs + opts(legend.position="none")
gs <- gs + labs(x=NULL,y="時間")
t = seq(7.5, 8, 1/6)
gs <- gs + scale_x_continuous(breaks=(sleep$年))
gs <- gs + scale_y_continuous(limit=c(7.5,8),breaks=t,labels=c("7:30", "7:40", "7:50", "8:00"))
gs <- gs + geom_text(aes(x=1996,y=value[年==1996],label=sex[年==1996]),vjust=-0.5)
gs <- gs + geom_text(aes(x=2006,y=value[年==2006],label=c("7:47", "7:32")),hjust=1,vjust=1.5)
gs
grid.gedit("text", gp=gpar(fontfamily="HiraKaku"))
パソコン世帯普及率
  • データ


    

  • グラフ
    • 

  • ソース
    • 






    
パソコン世帯普及率（内閣府2009年）

    

    

  • データ
    • 



    

  • グラフ
    • 

  • ソース
    • 






    
出生性比

    

    

  • データ
    • 



    

  • グラフ
    • 

  • ソース
    • 






    
出生率

    

    

  • データ
    • 



    

  • グラフ
    • 

  • ソース
    • 






    
生まれかわり

    

    

  • データ
    • 



    

  • グラフ
    • 

  • ソース
    • 






    
Twitterのフォロワー数

    

    

  • データ
    • 



    

  • グラフ
    • 

  • ソース
    • 






    
白血病による死亡数の推移

    

    

  • データ
    • 



    

  • グラフ
    • 

  • ソース
    • 








    
ストリップチャート・ボックスプロット

    

    

  • データ
    • 



    

  • グラフ
    • 

  • ソース
    • 




    
Michelson-Morley の実験

    

    

  • データ
    • 



    

  • グラフ
    • 

  • ソース
    • 








    
（Clevelandの）ドットプロット

    


    
都道府県別人口（ドットプロット）

    

    

  • データ
    • 



    

  • グラフ
    • 

  • ソース
    • 






    
ドットプロット

    

    

  • データ
    • 



    

  • グラフ
    • 

  • ソース
    • 








    
分位数プロット（quantile plot）

    


    
分位数プロット

    

    

  • データ
    • 



    

  • グラフ
    • 

  • ソース
    • 








    
散布図

    


    
都道府県別人口・面積

    

    

  • データ
    • 



    

  • グラフ
    • 

  • ソース
    • 






    
散布図（対話型）

    

    

  • データ
    • 



    

  • グラフ
    • 

  • ソース
    • 








    
その他

    


    
高校数学のグラフ

    

    

  • データ
    • 



    

  • グラフ
    • 

  • ソース
    • 








    
まとめ