T. Yoshiizumi - tjs_a_guide Diff
- Added parts are displayed like this.
- Deleted parts are displayed
like this.
単純集計に関する覚え書き／単一回答の集計
〜 rubyとRを用いた処理 〜
最終更新日:2013/05/292013/06/01
以下に掲げるドキュメントは、
[[tjs_a.zip|http://cup.sakura.ne.jp/tjs_a.zip]]
に同梱されている tjs_a.txt と同じ内容です。
＜はじめに＞
アンケート調査の集計に関わる機会があると、いつも、その場しのぎのスクリプトを書いて済ませてきました。そして、しばらくしてから再び別の調査の集計を行う際に、同じようにスクリプトを書くといった具合です。前の経験が蓄積されておらず、生かされないままでした。
統計的な分析より前に、単純な集計、数え上げでかなり苦労します。そのノウハウを何とか次に生かせるようにしたいと思いつつ、でも、次の機会になると前のやり方をみごとに忘れています。
ということで、なかなか整理がつかないままではありますが、少しノウハウらしきものを記しておこうと思います。なるべく Rの集計機能を利用するようにしてみます。
今回は、単一回答(どれか１つだけにまるをつける設問への回答)の集計を取り上げます。複数回答の集計や、人数・価格・長さなどの数値の集計は、別の機会に書きたいと思います。
ソフトウェアとしては R と ruby を用います。拙作 rrxwin(rubyからRを利用するためのライブラリ) も使います。
rrxwinは [[rrx100.zip|http://cup.sakura.ne.jp/rrx100.zip]] に一式含まれています。
サンプルのスクリプトを試した環境は次のとおり。
* Windows XP | VISTA
* ruby ver 1.8.7 | 1.9.3
* R ver 2.15.2
--------
{{toc_here}}
----
!1. 欠損値のない単一回答の集計
!!(1) 基本表の作成
○ 要点１：集計はtable関数。「table(dtf$sex, dtf$view)」
○ 要点２：集計表の項目順序の入れ換え「tbl[c("男","女"),c(……)]」
○ 要点３：集計時の項目タイトルの付加「table(性別=dtf$sex, 意見=dtf$view)」
回答者のID、sex(性別)、view(意見：賛成・中立・反対)が書かれたcsvファイル data01.csv があるとします。100人分のデータです。次のような内容です。
ID,sex,view
ID003,男,賛成
ID004,男,中立
ID005,女,反対
このcsvを素材にして、性別と意見のクロス集計表を作成します。
次のような表です。
|| ||賛成||中立||反対
||男||12||18||16
||女||24||18||12
上は回答者の人数を数え上げた基本表です。後から、これに総計欄を付加したり、構成比(パーセント)の表も作ったりします。
このような欠損値のない素材データの集計は容易です。Rで簡単に処理できます。
具体的には次のようなRプログラムです。tjs01.rbのRプログラム部分を抜粋します。
−−−− tjs01.rb 抜粋ここから
dtf <- read.csv("data01.csv", header=T)
tbl <- table(dtf$sex, dtf$view)
cat("基本表\n")
tbl
−−−− tjs01.rb 抜粋ここまで
Windows環境で ruby と R がインストールされている状態で、次のように実行します。
ruby.exe tjs01.rb [enter]
このように実行すると、基本表が標準出力に出力されます。
◇ 項目順序の入れ換え
tjs01.rbで作った基本表をみると、「女」が１行目、「男」が２行目になっています。この順序を入れ換える場合は次のようにします。
−−−− tjs01b.rb 抜粋ここから
dtf <- read.csv("data01.csv", header=T)
xx <- table(dtf$sex, dtf$view)
tbl <- xx[c("男", "女"),]
cat("基本表\n")
tbl
−−−− tjs01b.rb 抜粋ここまで
上のようにすると、「男」が１行目にきます。
もしviewの項目を反対・中立・賛成の順にしたければ、次のようにします。
tbl <- xx[c("男", "女"), c("反対", "中立", "賛成")]
◇ 項目タイトルの付加
今回のケースは、縦方向に「女・男」、横方向に「賛成・中立・反対」が列ぶので、項目タイトルを付加しなくても意味が分かります。しかし、項目タイトルがないと分かりにくいことがあります。
項目タイトルとして「性別」と「意見」を付加する例を掲げます。
−−−− tjs01c.rb 抜粋ここから
dtf <- read.csv("data01.csv", header=T)
tbl <- table(性別=dtf$sex, 意見=dtf$view)
cat("基本表\n")
tbl
−−−− tjs01b.rb 抜粋ここまで
上のRプログラムで、tblには項目タイトルの情報が付加されているわけですが、write.csv() などでそれを出力すると、項目タイトルは出力されません。
write.csv(tbl, file="output.csv")
上の記述で tbl の内容をoutput.csv に出力できますが、このcsvファイルには項目タイトルが付いてないので注意して下さい。
項目名と項目タイトルを表示したい場合は、dimnames(tbl) を用います。Rプログラムの中にこの１行を置くと、次のような表示が出ます。
−−−−
$性別
[1] "女" "男"
$意見
[1] "賛成" "中立" "反対"
−−−−
項目タイトルだけを取得したい時は names(dimnames(tbl)) です。この戻り値は２つの要素からなるベクトルで、第１要素が縦方向に流れる項目のタイトル、２番目が横方向に流れる項目のタイトルです。
参考まで、rubyの側でそれらの情報を取得するためのスクリプトを下に掲げておきます。
−−−− tjs01d.rb 抜粋ここから
rpro = <<EOS
dtf <- read.csv("data01.csv", header=T)
tbl <- table(性別=dtf$sex, 意見=dtf$view)
robj(tbl)
robj(names(dimnames(tbl)))
robj(dimnames(tbl))
EOS
##
hs = Rrx.rexec(rpro)
tbl = hs["robj1"]["self"]
print Rrx.ary2str(tbl, "\t")
printf("\n")
nn = hs["robj2"]["self"]
p nn
printf("\n")
xx =hs["robj2"]hs["robj3"]
xx.each do |key, val|
printf("%s\t", key)
p val
end
−−−− tjs01d.rb 抜粋ここまで
--------
!!(2) 総計欄の付加
○ 要点：総計欄の付加はaddmargins関数。「addmargins(tbl,1)」
Rプログラムでは addmargins() によって総計欄を付加できます。
基本表が tbl に入っている場合、次のように付加します。
addmargins(tbl, 1) # 行を付加(縦方向の総計欄を付加)
addmargins(tbl, 2) # 列を付加(横方向の総計欄を付加)
addmargins(tbl) # 行と列を付加(縦・横両方の総計欄を付加)
念のため Rプログラムとして掲げておきます。
−−−− tjs02.rb 抜粋ここから
dtf <- read.csv("data01.csv", header=T)
xx <- table(性別=dtf$sex, 意見=dtf$view)
tbl <- xx[c("男", "女"),]
cat("縦方向の総計欄を付加\n")
addmargins(tbl, 1)
cat("\n")
cat("横方向の総計欄を付加\n")
addmargins(tbl, 2)
cat("\n")
cat("縦・横両方の総計欄を付加\n")
addmargins(tbl)
−−−− tjs02.rb 抜粋ここまで
--------
!!(3) 構成比(パーセント)の算出
○ 要点：構成比はprop.table関数。「prop.table(tbl,1)*100」
Rプログラムにおいて、構成比は prop.table() で算出します。
基本表が tbl に入っている場合、次のようにします。
prop.table(tbl, 1)*100 # 横方向の総数に対する比率
prop.table(tbl, 2)*100 # 縦方向の総数に対する比率
構成比を算出する場合、総計欄を設ける必要がないのであれば、上の２つのパターンで十分です。
しかし、総計欄を付加したい時は少々面倒です。どのタイミングで総計欄を付加すべきかが分かりにくい。
文章で説明すると ややこしいですが、大雑把にいうと次のようになります。
まず、縦方向の総計欄(男子と女子の総数)を付加した上で、「男・女・Sum」の３つそれぞれについて横方向(賛成・中立・反対)の総数に対するパーセントを算出します。一番右端の縦の列は、どれも100%となります。
次に、横方向の総計欄(賛成・中立・反対の総数)を付加した上で、「賛成・中立・反対・Sum」の４つについて縦方向(男・女)の総数に対するパーセントを算出します。一番下端の横の列は、どれも100%です。
作られる表を示すと下のとおり。
横の総数に対する構成比の表
|| ||意見
||性別||賛成||中立||反対||Sum
||男||34.8||39.1||26.1||100
||女||22.2||33.3||44.4||100
||Sum||28.0||36.0||36.0||100
縦の総数に対する構成比の表
|| ||意見
||性別||賛成||中立||反対||Sum
||男||57.1||50||33.3||46
||女||42.9||50||66.7||54
||Sum||100.0||100||100.0||100
上の表を出力するためのRプログラムは次のとおり。
−−−− tjs03.rb 抜粋ここから
dtf <- read.csv("data01.csv", header=T)
xx <- table(性別=dtf$sex, 意見=dtf$view)
tbl <- xx[c("男", "女"),]
xx <- prop.table(addmargins(tbl,1), 1)*100
pct1 <- addmargins(xx, 2)
cat("横の総数に対する構成比の表\n")
round(pct1, 1)
cat("\n")
xx <- prop.table(addmargins(tbl,2), 2)*100
pct2 <- addmargins(xx, 1)
cat("縦の総数に対する構成比の表\n")
round(pct2, 1)
−−−− tjs03.rb 抜粋ここまで
--------
!!(4) 度数分布と構成比をrubyの側で得る
○ 要点：Rのlist関数で、複数のオブジェクトを一つにまとめる。
これまでのサンプルを取りまとめて、rubyの側で度数分布の表と構成比の表を得るrubyスクリプトを掲げておきます。拙作rrxwinを用います。
度数分布の表を１つ、構成比の表を２つ、総数３つの表を出力します。タブ区切り形式で、標準出力に出力します。
Windows環境で実行するとの前提で、rrxwinを用いていますが、linuxなどで実行する時は「rrxwin」を「rrx」に書き換えるなど修正して下さい。
−−−− tjs04.rb ここから
#! ruby -Ks
# list関数で複数のオブジェクトを一つにまとめる (coding: Windows-31J)
require "rrxwin"
##
rpro = <<EOS
dtf <- read.csv("data01.csv", header=T)
xx <- table(性別=dtf$sex, 意見=dtf$view)
tbl <- xx[c("男", "女"),]
xx <- prop.table(addmargins(tbl,1), 1)*100
pct1 <- addmargins(xx, 2)
xx <- prop.table(addmargins(tbl,2), 2)*100
pct2 <- addmargins(xx, 1)
robj(list(tbl=addmargins(tbl), pct1=round(pct1,1), pct2=round(pct2,1)))
EOS
##
hs = Rrx.rexec(rpro)
xx = hs["robj1"]
keys = %w(tbl pct1 pct2)
titles = %w(度数分布 横の総数に対する構成比 縦の総数に対する構成比)
keys.each do |key|
puts titles.shift
print Rrx.ary2str(xx[key], "\t")
printf("\n")
end
−−−− tjs04.rb ここまで
Rプログラムの中で、度数分布の表(tbl)と２つの構成比表(pct1, pct2)を１つのリストに収納し、robj() でruby側に渡せるようにしています。
ruby側では、「xx = hs["robj1"]」とした後で、xx["tbl"], xx["pct1"], xx["pct2"] で３つの表を参照できます。それぞれが２次元配列です。
--------
!!(5) 数字入力された素材データへの対応
○ 要点１：Rのベクトルなどの各要素の置換は「xx[xx==1] <- "男"」
○ 要点２：Rのfactorオブジェクトの二重性(整数値と文字列)に注意。
これまでの data01.csv は、性別欄に「男」 「女」と書かれていました。意見の欄には「賛成」 「中立」 「反対」が書かれています。
しかし、実際に入力作業する場合、これらの漢字を記入するよりも「男」の代わりに数字の 1、「女」の代わりに 2 を入力する方が楽です。同じように、賛成:1, 中立:2, 反対:3 の対応を想定して、半角の数字で入力すると楽です。
そのように数字に置き換えて入力してある素材データを仮に test.csv とします。その処理を考えます。
まず思いつくのは、素材データの数字を本来の漢字に置き換える方法です。ruby側で予め置き換えてからRプログラムで処理してもいいでしょうし、Rプログラムの read.csv() で素材データを読み込んでから置き換えを行ってもいいでしょう。
参考まで、後者の方法の一例を掲げてみます。test.csv を read.csv() で読み込んでから、置き換えを行った上で output.csv として書き出します。
−−−− Rプログラムここから
my.replace <- function(xx, elm) {
for (i in 1:length(elm)) {xx[xx==i] <- elm[i]}
return(xx)
}
dtf <- read.csv("test.csv", header=T)
dtf$sex <- my.replace(dtf$sex, c("男", "女"))
dtf$view <- my.replace(dtf$view, c("賛成", "中立", "反対"))
write.csv(dtf, file="output.csv")
−−−− Rプログラムここまで
Rプログラムでは、ベクトルなど複数の要素を持つオプジェクト xx について、その要素の置換に次のような記述を使えます。
xx[xx==1] <- "男" # 「1」をすべて「男」に置換
xx[xx==2] <- "女" # 「2」をすべて「女」に置換
上に掲げたRプログラムは、この置換方法を使っています。
置換のための記述を何度も書くのが面倒なので、my.replace() という関数を定義して対応しています。
さて、この方法には弱点があります。数字が示す順番の情報が無効になってしまう点です。
「男:1, 女:2」という対応にしていたわけですが、置換後の順序は逆になり、「女:1, 男:2」になります。
Rプログラムの内部処理で、文字コードの順に並べられるためそうなるのだと思います。
もう一つ弱点があります。前述した置換の方法は、Rプログラムの factor というオブジェクトに対し、素直には適用できないことです。
例えば、「男:1, 女:2」ではなく、「男:m, 女:f」というアルファベットでの代替を考えたとします。素材データの中に m, f があると、当然ながら数値でなく文字として認識されます。
では、その素材データを read.csv() で読み込んだ時に、dtf$sex が文字から構成されるvectorになるかというと、そうはなりません。factorになります。
「それがどうした？」と思われるかもしれませんが、factorであるが故に、次のような記述で置換することができません。
xx[xx=="m"] <- "男"
xx[xx=="f"] <- "女"
「男:1, 女:2」の場合は factor にならず、数値からなるvectorだったために上の置換方法を使えました。しかし、今度は駄目です。
この弱点を補いつつ置換を行うには、Rプログラムの factor の levels を調整します。
正確な言い方ではないかもしれませんが、factorというのは、Rの内部では整数値で管理されていて、各々の整数値に文字列が割り当てられている、そんなイメージです。全部を文字列で保持しているのに比べて、整数値で保持する方が記憶容量を節約できます。そして、おそらく処理速度の向上にもつながります。
いわば数値と文字列の二重性をおびているのが factor です。
xxがfactorである場合、levels(xx) は、文字列(整数値に対応する文字列)から構成されるベクトルを返します。整数値の小さい順番で列べられたものです。100人分のデータであっても、その欄に m, f の２種類しか書かれていなければ、levels(xx) の要素は２つだけです。
また、as.numeric(xx) は、整数値から構成されるベクトルを返します(というより変換します)。factorの要素と同じ個数だけ整数値が列べられたものになります。100人分のデータであれば 100個です。該当欄に m, f の２種類しか書かれていなくても(つまり整数値 1, 2 の２種類しかないとしても)、それが100人分列んだベクトルになります。
data02.csv は、性別欄に「男:m, 女:f」を記入し、意見欄に「賛成:1, 中立:2, 反対:3」を記入したデータです。
集計した時に、項目名がちゃんと漢字になるようにするには、mを男、1を賛成に書き換えたりするのではなく、数値と文字列の二重性をおびたfactorに変換します。既にfactorになっているものであっても、改めてfactorとして再設定します。
Rプログラムの中で factor() という関数を用いて、性別欄(dtf$sex)と意見欄(dtf$view)をこれにわたして調整します。次のように呼び出します。
dtf$sex <- factor(dtf$sex, levels=c("m","f"), labels=c("男","女"))
dtf$view <- factor(dtf$view, levels=1:3, labels=c("賛成","中立","反対"))
もし順番を「女・男」 「反対・中立・賛成」にしたければ次のようにします。
dtf$sex <- factor(dtf$sex, levels=c("f","m"), labels=c("女","男"))
dtf$view <- factor(dtf$view, levels=3:1, labels=c("反対","中立","賛成"))
なお、factor() にわたす引数の順番をかえないのであれば、「levels=」 「labels=」を省略してもいいようです。次のように書いて大丈夫です。
dtf$sex <- factor(dtf$sex, c("m","f"), c("男","女"))
dtf$view <- factor(dtf$view, 1:3, c("賛成","中立","反対"))
以下に、tjs04.rbと同じ出力結果をもたらすrubyスクリプト tjs05.rb を掲げます。素材データは data02.csv です。
度数分布と２つの構成比の表、総数３つの表を出力します。
−−−− tjs05.rb 抜粋ここから
rpro = <<EOS
dtf <- read.csv("data02.csv", header=T)
dtf$sex <- factor(dtf$sex, levels=c("m","f"), labels=c("男","女"))
dtf$view <- factor(dtf$view, levels=1:3, labels=c("賛成","中立","反対"))
tbl <- table(性別=dtf$sex, 意見=dtf$view)
xx <- prop.table(addmargins(tbl,1), 1)*100
pct1 <- addmargins(xx, 2)
xx <- prop.table(addmargins(tbl,2), 2)*100
pct2 <- addmargins(xx, 1)
robj(list(tbl=addmargins(tbl), pct1=round(pct1,1), pct2=round(pct2,1)))
EOS
##
hs = Rrx.rexec(rpro)
xx = hs["robj1"]
［後略］
−−−− tjs05.rb 抜粋ここまで
!!［参考］ 項目名の情報をファイルにメモしてある場合の処理
「男:m, 女:f」とか「賛成:1, 中立:2, 反対:3」といった対応は、メモしておかないと忘れます。
そこで、次のような内容のファイル data02.mem を作っておくことにします。２行からなるファイルです。
sex m,f 男,女
view 1,2,3 賛成,中立,反対
これは、data02.csvを見る時に、sexの欄では mが男、fが女を表し、viewの欄では 1が賛成、2が中立、3が反対を意味することを示します。
行の先頭がフィールド名、それ以降は、置換前の記号と置換後の名前をカンマ区切りで記します。
フィールド名、置換前の記号、置換後の名前の３つの区切りには半角スペースかタブコードを置きます。
このファイルを読み込んで、それを基にして tjs05.rb と同じ処理を行うことにします。
やり方はいろいろあると思いますが、rubyの側で、data02.memを読み込んで、Rプログラムの一部を書き出すという方法を採ります。すべてRで行うのも可能だと思いますが、私には少々敷居が高いので rubyを用います。
rubyの側で remake_factor() というメソッドを定義します。これは、次のように用います。
remake_factor("data02.mem", "dtf")
第１引数がファイル名、第２引数は Rプログラムにおけるデータフレームの名前です。
このメソッドの戻り値は下の２行です。
dtf$sex <- factor(dtf$sex, levels=c("m","f"), labels=c("男","女"))
{{br}}
dtf$view <- factor(dtf$view, levels=c("1","2","3"), labels=c("賛成","中立","反対"))
この２行をRプログラムの一部として挿入します。tjs05b.rb では erb を用いて挿入しています。
tjs05b.rb でやっていることは簡単なのですが、ちょっとごちゃごちゃしています。少しでも整理するため、rubyのメソッド定義を tjsf.rb に書き込んでおいて、それを取り込む形にしてあります。
erbを利用するサンプルにもなっているので、少し長いですが下に tjs05b.rb を全部掲げます。
−−−− tjs05b.rb ここから
#! ruby -Ks
# メモファイルに基づく数値・記号の置換＆erbの利用 (coding: Windows-31J)
require "erb"
require "rrxwin"
require "./tjsf"
##
rpro = <<EOS
dtf <- read.csv("<%= csv_file %>", header=T)
<%= remake_factor(memo_file, "dtf") %>
tbl <- table(性別=dtf$sex, 意見=dtf$view)
xx <- prop.table(addmargins(tbl,1), 1)*100
pct1 <- addmargins(xx, 2)
xx <- prop.table(addmargins(tbl,2), 2)*100
pct2 <- addmargins(xx, 1)
robj(list(tbl=addmargins(tbl), pct1=round(pct1,1), pct2=round(pct2,1)))
EOS
##
csv_file = "data02.csv"
memo_file = "data02.mem"
erb = ERB.new(rpro)
hs = Rrx.rexec(erb.result(binding))
xx = hs["robj1"]
keys = %w(tbl pct1 pct2)
titles = %w(度数分布 横の総数に対する構成比 縦の総数に対する構成比)
keys.each do |key|
puts titles.shift
print Rrx.ary2str(xx[key], "\t")
printf("\n")
end
−−−− tjs05b.rb ここまで
--------
!2. 欠損値のある単一回答の集計
回答者のID、sex(性別)、view(意見：賛成・中立・反対)が書かれたcsvファイル data03.csv を素材にします。
data02.csvと概ね同じですが、性別または意見の欄が空欄のケースがあります。いわば無回答があるわけです。
無回答がある場合の集計方法として、1)有効回答だけに着目する方法(無回答を除いた有効回答総数を100%とみなす)、2)無回答も１種の回答とみなす方法(無回答を含めた総数を100%とする)の２通りがあります。
どちらを採用するかは、調査項目の内容や、「無回答の背景にあると思われるもの」の判断の仕方によって違ってきます。
ということで、２通りの集計方法それぞれを取り上げてみます。
!!(1) 無回答を含めた集計
○ 要点１：欠損値を含む集計は「table(dtf$sex, dtf$view, exclude=NULL)」
○ 要点２：read.csv関数における欠損値の指定は na.strings オプション
有効回答のみに着目した集計、つまり無回答に全く言及しない集計は簡単です。tjs05.rbをそのまま使えます。
集計した度数分布は次のようになります。
|| ||賛成||中立||反対||Sum
||男||14||14||11||39
||女||9||17||20||46
||Sum||23||31||31||85
無回答が僅かで、考慮しなくてもいい旨を文章で説明していれば、上の表でも大丈夫かもしれません。
しかし、一般的には情報不足です。文章を読まずに表だけ拾い読みすることがあるので、やはり無回答がどれくらいかを表に盛り込む方が望ましいといえます。
Rプログラムでは、欠損値をNAという特別の値として扱うことが多いようです。NA は Not Available の意味。
table() 関数は、通常、NAを除いて集計しますが、「exclude=NULL」というオプションを指定すると、NAを含めた集計を行います。例えば次のようにします。
tbl <- table(性別=dtf$sex, 意見=dtf$view, exclude=NULL)
data03.csvの集計結果(度数分布)は下のとおりです。
|| ||賛成||中立||反対||NA||Sum
||男||14||14||11||4||43
||女||9||17||20||4||50
||NA||2||3||2||0||7
||Sum||25||34||33||8||100
「NA」の項をみると、無回答がどのくらいあるか把握できます。
以下にスクリプト tjs06.rb を示します。
−−−− tjs06.rb 抜粋ここから
rpro = <<EOS
dtf <- read.csv("data03.csv", header=T, na.strings="")
dtf$sex <- factor(dtf$sex, levels=c("m","f"), labels=c("男","女"))
dtf$view <- factor(dtf$view, levels=1:3, labels=c("賛成","中立","反対"))
tbl <- table(性別=dtf$sex, 意見=dtf$view, exclude=NULL)
xx <- prop.table(addmargins(tbl,1), 1)*100
pct1 <- addmargins(xx, 2)
xx <- prop.table(addmargins(tbl,2), 2)*100
pct2 <- addmargins(xx, 1)
robj(list(tbl=addmargins(tbl), pct1=round(pct1,1), pct2=round(pct2,1)))
EOS
##
hs = Rrx.rexec(rpro)
xx = hs["robj1"]
［後略］
−−−− tjs06.rb 抜粋ここまで
Rプログラムの read.csv() の最後の引数として「na.strings=""」を記述しています。これは、csvの空欄をNAにするとの指定です。
実は、この記述はなくてもいいのですが、念のため明示してみました。
もしcsvデータの中に NONE という４文字を書き込んでおいて、それを無効なものとして扱うのであれば、「na.strings="NONE"」とします。
性別欄は「m, f」の文字なので NONE にするとしても、意見の欄は「1, 2, 3」なので -999 をNAとして書き込みたい、ということであれば、次のようにします。
read.csv("data03.csv", header=T, na.strings=c("NONE",-999))
この場合、NONEがどこにあっても(性別欄以外のところであっても) NAになります。同じように、-999がどこにあっても(意見欄以外であっても) NAになるので注意して下さい。
!!［補足］ 行と列の項目名を変更する(「NA, Sum」→「無回答, 総数」)
tjs06.rb で作成した表には、NA, Sum という項目名がありました。これを「無回答」 「総数」に変更します。
xxが行列(tableやmatrix)である場合、colnames(xx) で列の項目名、rownames(xx) で行の項目名を得られます。これら両関数の戻り値は vector で、項目名が割り当てられていなければ NULL を返します。
列の項目名の中の NA を「無回答」に変更するには次のようにします。
colnames(xx)[is.na(colnames(xx))] <- "無回答"
また、列の項目名の中の「Sum」を「総数」に変更する場合は次のとおり。
colnames(xx)[colnames(xx)=="Sum"] <- "総数"
上の colnames を rownames にすれば、行の項目名を変更することになります。
これら変更の手続きをRプログラムの中に逐一 書き込むのは面倒なので、remake.item() という関数定義を tjsf.r に盛り込んで、それを呼び出すことにします。
xxが行列である場合、次のように用います。いくつかのパターンを示します。
xx <- remake.item(xx, c(NA,"Sum"), c("無回答","総数"))
xx <- remake.item(xx) # 「NA, Sum」→「無回答, 総数」
xx <- remake.item(xx, c("無効回答","総計")) # 「無効回答, 総計」に変更
基本的には３つの引数を指定します。第２引数は旧い名前群、第３引数は新しい名前群です。行と列の両方について、旧い名前群を新しい名前群に置き換えます。
第３引数を省略すると、第２引数が新しい名前群とみなされ、旧い名前群はデフォルト値の c(NA, "Sum") であるものとされます。
第２引数と第３引数の両方を省略すると、それぞれデフォルト値の c(NA, "Sum"), c("無回答", "総数") が指定されたものとみなされます。
--------
!!(2) 行列から行・列を取り出して引き算する(有効回答数の算出)
○ 要点：行列の5列目から4列目を引き算するには「xx[,5] - xx[,4]」
無回答を含む基本表(度数分布)を改めて示します。基本的に、tjs06.rbで作成される表と同じ内容です。
|| ||賛成||中立||反対||無回答||総数
||男||14||14||11||4||43
||女||9||17||20||4||50
||無回答||2||3||2||0||7
||総数||25||34||33||8||100
上の表をみても、有効回答がどれくらいかすぐには分かりません。
「総数」から「無回答」を引き算してやれば、有効回答を算出できるのでそれを試みます。とりあえず、縦方向に流れる列に着目します。
行列から特定の列を取り出すには xx[,4] とか xx[,5] とします。「無回答」が４列目、「総数」が５列目なので、その番号で取り出します。
また、項目名で取り出すこともできます。xx[,"無回答"], xx[,"総数"] と書くことができます。
つまり、有効回答は次のようにして求めることができます。２通りを示します。
yy <- xx[,5] - xx[,4]
yy <- xx[,"総数"] - xx[,"無回答"]
一方、縦方向の列ではなく、横方向の行を取り出して引き算する場合は次のとおりです。
行に着目すると、「無回答」が３番目、「総数」が４番目です。
yy <- xx[4,] - xx[3,]
yy <- xx["総数",] - xx["無回答",]
ここまで述べてきたノウハウを用いて、有効回答と無回答の両方を含む度数分布の表を作成してみます。
行にしろ列にしろ、「無回答」と「総数」が最後の２つの項目(つまり外側に位置する２つ)であることを利用して、行と列の取り出しを番号で行います。
また、行列の結合には cbind(), rbind() を使います。
以下、tjs07.rb です。
−−−− tjs07.rb 抜粋ここから
rpro = <<EOS
source("tjsf.r")
dtf <- read.csv("data03.csv", header=T, na.strings="")
dtf$sex <- factor(dtf$sex, levels=c("m","f"), labels=c("男","女"))
dtf$view <- factor(dtf$view, levels=1:3, labels=c("賛成","中立","反対"))
tbl <- table(性別=dtf$sex, 意見=dtf$view, exclude=NULL)
xx <- addmargins(tbl) # 行・列の両方に総数を付加
xx <- remake.item(xx) # 「NA, Sum」→「無回答, 総数」
n <- ncol(xx) # 列の個数
yy <- xx[,n] - xx[,n-1] # 有効回答数を算出
zz <- cbind(xx[,1:(n-2)], 有効回答=yy, xx[,(n-1):n])
xx <- zz
n <- nrow(xx) # 行の個数
yy <- xx[n,] - xx[n-1,] # 有効回答数を算出
zz <- rbind(xx[1:(n-2),], 有効回答=yy, xx[(n-1):n,])
robj(zz)
EOS
##
hs = Rrx.rexec(rpro)
zz = hs["robj1"]["self"]
print Rrx.ary2str(zz, "\t")
−−−− tjs07.rb 抜粋ここまで
cbind(), rbind() の第２引数を「有効回答=yy」としています。これを単に「yy」にすると、作成される表の項目名が「有効回答」にならず「yy」になります。
作成される表は次のとおりです。
|| ||賛成||中立||反対||有効回答||無回答||総数
||男||14||14||11||39||4||43
||女||9||17||20||46||4||50
||有効回答||23||31||31||85||8||93
||無回答||2||3||2||7||0||7
||総数||25||34||33||92||8||100
これで、有効回答がどれくらいか分かるようになりました。縦の有効回答と横の有効回答が交差する欄は 85 ですが、これは、性別と意見の両方とも分かる人の数です。
有効回答と無回答が交差する欄が何を意味するのか、ちょっと分かりにくい感じがありますが、性別が分かっているけれども意見が分からない人の数、その逆に、性別が分からないけれども意見が分かる人の数を示しています。
横の無回答と縦の無回答が交差する欄は、性別も意見も両方とも分からない人の数を示しますが、表ではそれが 0 になっています。
tjs07.rbのRプログラムの部分については、コメントを見ていただければ、どんなふうに処理を進めているのか分かると思います。
詳しい仕様は省略しますが、cbind() は、行列に列を追加する関数です。a, b が行列やベクトルである場合、「x<-cbind(a,b)」とすれば、aにbを追加した結果がxに入ります。xは、列を追加した結果なので、aよりも横に拡がったものになります。cbind(a,b,c) のように３つ以上の引数を与えることもできます。
rbind() は、行を追加するものです。追加した結果は、当然ながら縦に拡がります。
--------
!!(3) 総数に対する比率(パーセント)の表を作成する
○ 要点：１行目について最右端の列で割り算した結果(パーセント)を得るには{{br}}
「xx[1,] / xx[1,ncol(xx)] * 100」
tjs07.rbで作成した表は、度数分布の表です。これを基にして、各欄の値が「総数」に占める比率(パーセント)を求めてみます。そうすれば、有効回答率等が分かります。
ちなみに、有効回答を100%とする比率(構成比)は、tjs05.rbが出力する表で把握できます。
構成比を求めるのに、前のサンプルでは prop.table() を使いました。tjs03.rb がそれです。
しかし、今回は prop.table を利用しません。「有効回答」が一種の小計になっているため、「総数」が各欄を足し算したものにならないからです。
例えば、「男」に着目して、「賛成・中立・反対・有効回答・無回答」を素直に足し算すると、82 になります。これは「総数」の 43 よりずっと大きい数です。「有効回答」の分だけ多くなっています。
そこで、各欄を「総数」で割り算して比率を求めることにします。
とりあえず第１行目の「男」に着目すると、「総数」が最右端の６列目なので、「総数」を含めた各欄の比率(パーセント)は次のようにして算出できます。
male <- xx[1,] / xx[1,6] * 100
こうすると、「賛成・中立・反対・有効回答・無回答・総数」の各パーセントの値が male に入ります。
「男」だけでなく「女・有効回答・無回答・総数」の５行分について比率を求めようとすれば、次のようになります。基本表は xx に入っているものとします。
nc <- ncol(xx) # 最右端の列の番号
for (i in 1:nrow(xx)) {
if (xx[i,nc] == 0) next # 0による割り算は避ける
xx[i,] <- xx[i,] / xx[i,nc] * 100 # i行目を比率に入れ換え
}
上の処理を経ると、度数が記録されていたxxの中身が、比率(パーセント)に入れ換えられます。
tjs07.rbで得られる基本表に対し、上の処理を施して比率を得るようにしたのが tjs08.rb です。
参考まで掲げておきます。
−−−− tjs08.rb 抜粋ここから
rpro = <<EOS
source("tjsf.r")
dtf <- read.csv("data03.csv", header=T, na.strings="")
dtf$sex <- factor(dtf$sex, levels=c("m","f"), labels=c("男","女"))
dtf$view <- factor(dtf$view, levels=1:3, labels=c("賛成","中立","反対"))
tbl <- table(性別=dtf$sex, 意見=dtf$view, exclude=NULL)
xx <- addmargins(tbl) # 行・列の両方に総数を付加
xx <- remake.item(xx) # 「NA, Sum」→「無回答, 総数」
n <- ncol(xx) # 列の個数
yy <- xx[,n] - xx[,n-1] # 有効回答数を算出
zz <- cbind(xx[,1:(n-2)], 有効回答=yy, xx[,(n-1):n])
xx <- zz
n <- nrow(xx) # 行の個数
yy <- xx[n,] - xx[n-1,] # 有効回答数を算出
tbl <- rbind(xx[1:(n-2),], 有効回答=yy, xx[(n-1):n,])
nc <- ncol(tbl)
nr <- nrow(tbl)
pct1 <- tbl
for (i in 1:nr) {
if (pct1[i,nc] == 0) next
pct1[i,] <- pct1[i,] / pct1[i,nc] * 100
}
pct2 <- tbl
for (i in 1:nc) {
if (pct2[nr,i] == 0) next
pct2[,i] <- pct2[,i] / pct2[nr,i] * 100
}
robj(list(tbl=tbl, pct1=round(pct1,1), pct2=round(pct2,1)))
EOS
##
hs = Rrx.rexec(rpro)
xx = hs["robj1"]
［後略］
−−−− tjs08.rb 抜粋ここまで
--------
以上で「単純集計に関する覚え書き／単一回答の集計」をおわりにします。
Copyright (C) T. Yoshiizumi, 2013 All rights reserved.
〜 rubyとRを用いた処理 〜
最終更新日:
以下に掲げるドキュメントは、
[[tjs_a.zip|http://cup.sakura.ne.jp/tjs_a.zip]]
に同梱されている tjs_a.txt と同じ内容です。
＜はじめに＞
アンケート調査の集計に関わる機会があると、いつも、その場しのぎのスクリプトを書いて済ませてきました。そして、しばらくしてから再び別の調査の集計を行う際に、同じようにスクリプトを書くといった具合です。前の経験が蓄積されておらず、生かされないままでした。
統計的な分析より前に、単純な集計、数え上げでかなり苦労します。そのノウハウを何とか次に生かせるようにしたいと思いつつ、でも、次の機会になると前のやり方をみごとに忘れています。
ということで、なかなか整理がつかないままではありますが、少しノウハウらしきものを記しておこうと思います。なるべく Rの集計機能を利用するようにしてみます。
今回は、単一回答(どれか１つだけにまるをつける設問への回答)の集計を取り上げます。複数回答の集計や、人数・価格・長さなどの数値の集計は、別の機会に書きたいと思います。
ソフトウェアとしては R と ruby を用います。拙作 rrxwin(rubyからRを利用するためのライブラリ) も使います。
rrxwinは [[rrx100.zip|http://cup.sakura.ne.jp/rrx100.zip]] に一式含まれています。
サンプルのスクリプトを試した環境は次のとおり。
* Windows XP | VISTA
* ruby ver 1.8.7 | 1.9.3
* R ver 2.15.2
--------
{{toc_here}}
----
!1. 欠損値のない単一回答の集計
!!(1) 基本表の作成
○ 要点１：集計はtable関数。「table(dtf$sex, dtf$view)」
○ 要点２：集計表の項目順序の入れ換え「tbl[c("男","女"),c(……)]」
○ 要点３：集計時の項目タイトルの付加「table(性別=dtf$sex, 意見=dtf$view)」
回答者のID、sex(性別)、view(意見：賛成・中立・反対)が書かれたcsvファイル data01.csv があるとします。100人分のデータです。次のような内容です。
ID,sex,view
ID003,男,賛成
ID004,男,中立
ID005,女,反対
このcsvを素材にして、性別と意見のクロス集計表を作成します。
次のような表です。
|| ||賛成||中立||反対
||男||12||18||16
||女||24||18||12
上は回答者の人数を数え上げた基本表です。後から、これに総計欄を付加したり、構成比(パーセント)の表も作ったりします。
このような欠損値のない素材データの集計は容易です。Rで簡単に処理できます。
具体的には次のようなRプログラムです。tjs01.rbのRプログラム部分を抜粋します。
−−−− tjs01.rb 抜粋ここから
dtf <- read.csv("data01.csv", header=T)
tbl <- table(dtf$sex, dtf$view)
cat("基本表\n")
tbl
−−−− tjs01.rb 抜粋ここまで
Windows環境で ruby と R がインストールされている状態で、次のように実行します。
ruby.exe tjs01.rb [enter]
このように実行すると、基本表が標準出力に出力されます。
◇ 項目順序の入れ換え
tjs01.rbで作った基本表をみると、「女」が１行目、「男」が２行目になっています。この順序を入れ換える場合は次のようにします。
−−−− tjs01b.rb 抜粋ここから
dtf <- read.csv("data01.csv", header=T)
xx <- table(dtf$sex, dtf$view)
tbl <- xx[c("男", "女"),]
cat("基本表\n")
tbl
−−−− tjs01b.rb 抜粋ここまで
上のようにすると、「男」が１行目にきます。
もしviewの項目を反対・中立・賛成の順にしたければ、次のようにします。
tbl <- xx[c("男", "女"), c("反対", "中立", "賛成")]
◇ 項目タイトルの付加
今回のケースは、縦方向に「女・男」、横方向に「賛成・中立・反対」が列ぶので、項目タイトルを付加しなくても意味が分かります。しかし、項目タイトルがないと分かりにくいことがあります。
項目タイトルとして「性別」と「意見」を付加する例を掲げます。
−−−− tjs01c.rb 抜粋ここから
dtf <- read.csv("data01.csv", header=T)
tbl <- table(性別=dtf$sex, 意見=dtf$view)
cat("基本表\n")
tbl
−−−− tjs01b.rb 抜粋ここまで
上のRプログラムで、tblには項目タイトルの情報が付加されているわけですが、write.csv() などでそれを出力すると、項目タイトルは出力されません。
write.csv(tbl, file="output.csv")
上の記述で tbl の内容をoutput.csv に出力できますが、このcsvファイルには項目タイトルが付いてないので注意して下さい。
項目名と項目タイトルを表示したい場合は、dimnames(tbl) を用います。Rプログラムの中にこの１行を置くと、次のような表示が出ます。
−−−−
$性別
[1] "女" "男"
$意見
[1] "賛成" "中立" "反対"
−−−−
項目タイトルだけを取得したい時は names(dimnames(tbl)) です。この戻り値は２つの要素からなるベクトルで、第１要素が縦方向に流れる項目のタイトル、２番目が横方向に流れる項目のタイトルです。
参考まで、rubyの側でそれらの情報を取得するためのスクリプトを下に掲げておきます。
−−−− tjs01d.rb 抜粋ここから
rpro = <<EOS
dtf <- read.csv("data01.csv", header=T)
tbl <- table(性別=dtf$sex, 意見=dtf$view)
robj(tbl)
robj(names(dimnames(tbl)))
robj(dimnames(tbl))
EOS
##
hs = Rrx.rexec(rpro)
tbl = hs["robj1"]["self"]
print Rrx.ary2str(tbl, "\t")
printf("\n")
nn = hs["robj2"]["self"]
p nn
printf("\n")
xx =
xx.each do |key, val|
printf("%s\t", key)
p val
end
−−−− tjs01d.rb 抜粋ここまで
--------
!!(2) 総計欄の付加
○ 要点：総計欄の付加はaddmargins関数。「addmargins(tbl,1)」
Rプログラムでは addmargins() によって総計欄を付加できます。
基本表が tbl に入っている場合、次のように付加します。
addmargins(tbl, 1) # 行を付加(縦方向の総計欄を付加)
addmargins(tbl, 2) # 列を付加(横方向の総計欄を付加)
addmargins(tbl) # 行と列を付加(縦・横両方の総計欄を付加)
念のため Rプログラムとして掲げておきます。
−−−− tjs02.rb 抜粋ここから
dtf <- read.csv("data01.csv", header=T)
xx <- table(性別=dtf$sex, 意見=dtf$view)
tbl <- xx[c("男", "女"),]
cat("縦方向の総計欄を付加\n")
addmargins(tbl, 1)
cat("\n")
cat("横方向の総計欄を付加\n")
addmargins(tbl, 2)
cat("\n")
cat("縦・横両方の総計欄を付加\n")
addmargins(tbl)
−−−− tjs02.rb 抜粋ここまで
--------
!!(3) 構成比(パーセント)の算出
○ 要点：構成比はprop.table関数。「prop.table(tbl,1)*100」
Rプログラムにおいて、構成比は prop.table() で算出します。
基本表が tbl に入っている場合、次のようにします。
prop.table(tbl, 1)*100 # 横方向の総数に対する比率
prop.table(tbl, 2)*100 # 縦方向の総数に対する比率
構成比を算出する場合、総計欄を設ける必要がないのであれば、上の２つのパターンで十分です。
しかし、総計欄を付加したい時は少々面倒です。どのタイミングで総計欄を付加すべきかが分かりにくい。
文章で説明すると ややこしいですが、大雑把にいうと次のようになります。
まず、縦方向の総計欄(男子と女子の総数)を付加した上で、「男・女・Sum」の３つそれぞれについて横方向(賛成・中立・反対)の総数に対するパーセントを算出します。一番右端の縦の列は、どれも100%となります。
次に、横方向の総計欄(賛成・中立・反対の総数)を付加した上で、「賛成・中立・反対・Sum」の４つについて縦方向(男・女)の総数に対するパーセントを算出します。一番下端の横の列は、どれも100%です。
作られる表を示すと下のとおり。
横の総数に対する構成比の表
|| ||意見
||性別||賛成||中立||反対||Sum
||男||34.8||39.1||26.1||100
||女||22.2||33.3||44.4||100
||Sum||28.0||36.0||36.0||100
縦の総数に対する構成比の表
|| ||意見
||性別||賛成||中立||反対||Sum
||男||57.1||50||33.3||46
||女||42.9||50||66.7||54
||Sum||100.0||100||100.0||100
上の表を出力するためのRプログラムは次のとおり。
−−−− tjs03.rb 抜粋ここから
dtf <- read.csv("data01.csv", header=T)
xx <- table(性別=dtf$sex, 意見=dtf$view)
tbl <- xx[c("男", "女"),]
xx <- prop.table(addmargins(tbl,1), 1)*100
pct1 <- addmargins(xx, 2)
cat("横の総数に対する構成比の表\n")
round(pct1, 1)
cat("\n")
xx <- prop.table(addmargins(tbl,2), 2)*100
pct2 <- addmargins(xx, 1)
cat("縦の総数に対する構成比の表\n")
round(pct2, 1)
−−−− tjs03.rb 抜粋ここまで
--------
!!(4) 度数分布と構成比をrubyの側で得る
○ 要点：Rのlist関数で、複数のオブジェクトを一つにまとめる。
これまでのサンプルを取りまとめて、rubyの側で度数分布の表と構成比の表を得るrubyスクリプトを掲げておきます。拙作rrxwinを用います。
度数分布の表を１つ、構成比の表を２つ、総数３つの表を出力します。タブ区切り形式で、標準出力に出力します。
Windows環境で実行するとの前提で、rrxwinを用いていますが、linuxなどで実行する時は「rrxwin」を「rrx」に書き換えるなど修正して下さい。
−−−− tjs04.rb ここから
#! ruby -Ks
# list関数で複数のオブジェクトを一つにまとめる (coding: Windows-31J)
require "rrxwin"
##
rpro = <<EOS
dtf <- read.csv("data01.csv", header=T)
xx <- table(性別=dtf$sex, 意見=dtf$view)
tbl <- xx[c("男", "女"),]
xx <- prop.table(addmargins(tbl,1), 1)*100
pct1 <- addmargins(xx, 2)
xx <- prop.table(addmargins(tbl,2), 2)*100
pct2 <- addmargins(xx, 1)
robj(list(tbl=addmargins(tbl), pct1=round(pct1,1), pct2=round(pct2,1)))
EOS
##
hs = Rrx.rexec(rpro)
xx = hs["robj1"]
keys = %w(tbl pct1 pct2)
titles = %w(度数分布 横の総数に対する構成比 縦の総数に対する構成比)
keys.each do |key|
puts titles.shift
print Rrx.ary2str(xx[key], "\t")
printf("\n")
end
−−−− tjs04.rb ここまで
Rプログラムの中で、度数分布の表(tbl)と２つの構成比表(pct1, pct2)を１つのリストに収納し、robj() でruby側に渡せるようにしています。
ruby側では、「xx = hs["robj1"]」とした後で、xx["tbl"], xx["pct1"], xx["pct2"] で３つの表を参照できます。それぞれが２次元配列です。
--------
!!(5) 数字入力された素材データへの対応
○ 要点１：Rのベクトルなどの各要素の置換は「xx[xx==1] <- "男"」
○ 要点２：Rのfactorオブジェクトの二重性(整数値と文字列)に注意。
これまでの data01.csv は、性別欄に「男」 「女」と書かれていました。意見の欄には「賛成」 「中立」 「反対」が書かれています。
しかし、実際に入力作業する場合、これらの漢字を記入するよりも「男」の代わりに数字の 1、「女」の代わりに 2 を入力する方が楽です。同じように、賛成:1, 中立:2, 反対:3 の対応を想定して、半角の数字で入力すると楽です。
そのように数字に置き換えて入力してある素材データを仮に test.csv とします。その処理を考えます。
まず思いつくのは、素材データの数字を本来の漢字に置き換える方法です。ruby側で予め置き換えてからRプログラムで処理してもいいでしょうし、Rプログラムの read.csv() で素材データを読み込んでから置き換えを行ってもいいでしょう。
参考まで、後者の方法の一例を掲げてみます。test.csv を read.csv() で読み込んでから、置き換えを行った上で output.csv として書き出します。
−−−− Rプログラムここから
my.replace <- function(xx, elm) {
for (i in 1:length(elm)) {xx[xx==i] <- elm[i]}
return(xx)
}
dtf <- read.csv("test.csv", header=T)
dtf$sex <- my.replace(dtf$sex, c("男", "女"))
dtf$view <- my.replace(dtf$view, c("賛成", "中立", "反対"))
write.csv(dtf, file="output.csv")
−−−− Rプログラムここまで
Rプログラムでは、ベクトルなど複数の要素を持つオプジェクト xx について、その要素の置換に次のような記述を使えます。
xx[xx==1] <- "男" # 「1」をすべて「男」に置換
xx[xx==2] <- "女" # 「2」をすべて「女」に置換
上に掲げたRプログラムは、この置換方法を使っています。
置換のための記述を何度も書くのが面倒なので、my.replace() という関数を定義して対応しています。
さて、この方法には弱点があります。数字が示す順番の情報が無効になってしまう点です。
「男:1, 女:2」という対応にしていたわけですが、置換後の順序は逆になり、「女:1, 男:2」になります。
Rプログラムの内部処理で、文字コードの順に並べられるためそうなるのだと思います。
もう一つ弱点があります。前述した置換の方法は、Rプログラムの factor というオブジェクトに対し、素直には適用できないことです。
例えば、「男:1, 女:2」ではなく、「男:m, 女:f」というアルファベットでの代替を考えたとします。素材データの中に m, f があると、当然ながら数値でなく文字として認識されます。
では、その素材データを read.csv() で読み込んだ時に、dtf$sex が文字から構成されるvectorになるかというと、そうはなりません。factorになります。
「それがどうした？」と思われるかもしれませんが、factorであるが故に、次のような記述で置換することができません。
xx[xx=="m"] <- "男"
xx[xx=="f"] <- "女"
「男:1, 女:2」の場合は factor にならず、数値からなるvectorだったために上の置換方法を使えました。しかし、今度は駄目です。
この弱点を補いつつ置換を行うには、Rプログラムの factor の levels を調整します。
正確な言い方ではないかもしれませんが、factorというのは、Rの内部では整数値で管理されていて、各々の整数値に文字列が割り当てられている、そんなイメージです。全部を文字列で保持しているのに比べて、整数値で保持する方が記憶容量を節約できます。そして、おそらく処理速度の向上にもつながります。
いわば数値と文字列の二重性をおびているのが factor です。
xxがfactorである場合、levels(xx) は、文字列(整数値に対応する文字列)から構成されるベクトルを返します。整数値の小さい順番で列べられたものです。100人分のデータであっても、その欄に m, f の２種類しか書かれていなければ、levels(xx) の要素は２つだけです。
また、as.numeric(xx) は、整数値から構成されるベクトルを返します(というより変換します)。factorの要素と同じ個数だけ整数値が列べられたものになります。100人分のデータであれば 100個です。該当欄に m, f の２種類しか書かれていなくても(つまり整数値 1, 2 の２種類しかないとしても)、それが100人分列んだベクトルになります。
data02.csv は、性別欄に「男:m, 女:f」を記入し、意見欄に「賛成:1, 中立:2, 反対:3」を記入したデータです。
集計した時に、項目名がちゃんと漢字になるようにするには、mを男、1を賛成に書き換えたりするのではなく、数値と文字列の二重性をおびたfactorに変換します。既にfactorになっているものであっても、改めてfactorとして再設定します。
Rプログラムの中で factor() という関数を用いて、性別欄(dtf$sex)と意見欄(dtf$view)をこれにわたして調整します。次のように呼び出します。
dtf$sex <- factor(dtf$sex, levels=c("m","f"), labels=c("男","女"))
dtf$view <- factor(dtf$view, levels=1:3, labels=c("賛成","中立","反対"))
もし順番を「女・男」 「反対・中立・賛成」にしたければ次のようにします。
dtf$sex <- factor(dtf$sex, levels=c("f","m"), labels=c("女","男"))
dtf$view <- factor(dtf$view, levels=3:1, labels=c("反対","中立","賛成"))
なお、factor() にわたす引数の順番をかえないのであれば、「levels=」 「labels=」を省略してもいいようです。次のように書いて大丈夫です。
dtf$sex <- factor(dtf$sex, c("m","f"), c("男","女"))
dtf$view <- factor(dtf$view, 1:3, c("賛成","中立","反対"))
以下に、tjs04.rbと同じ出力結果をもたらすrubyスクリプト tjs05.rb を掲げます。素材データは data02.csv です。
度数分布と２つの構成比の表、総数３つの表を出力します。
−−−− tjs05.rb 抜粋ここから
rpro = <<EOS
dtf <- read.csv("data02.csv", header=T)
dtf$sex <- factor(dtf$sex, levels=c("m","f"), labels=c("男","女"))
dtf$view <- factor(dtf$view, levels=1:3, labels=c("賛成","中立","反対"))
tbl <- table(性別=dtf$sex, 意見=dtf$view)
xx <- prop.table(addmargins(tbl,1), 1)*100
pct1 <- addmargins(xx, 2)
xx <- prop.table(addmargins(tbl,2), 2)*100
pct2 <- addmargins(xx, 1)
robj(list(tbl=addmargins(tbl), pct1=round(pct1,1), pct2=round(pct2,1)))
EOS
##
hs = Rrx.rexec(rpro)
xx = hs["robj1"]
［後略］
−−−− tjs05.rb 抜粋ここまで
!!［参考］ 項目名の情報をファイルにメモしてある場合の処理
「男:m, 女:f」とか「賛成:1, 中立:2, 反対:3」といった対応は、メモしておかないと忘れます。
そこで、次のような内容のファイル data02.mem を作っておくことにします。２行からなるファイルです。
sex m,f 男,女
view 1,2,3 賛成,中立,反対
これは、data02.csvを見る時に、sexの欄では mが男、fが女を表し、viewの欄では 1が賛成、2が中立、3が反対を意味することを示します。
行の先頭がフィールド名、それ以降は、置換前の記号と置換後の名前をカンマ区切りで記します。
フィールド名、置換前の記号、置換後の名前の３つの区切りには半角スペースかタブコードを置きます。
このファイルを読み込んで、それを基にして tjs05.rb と同じ処理を行うことにします。
やり方はいろいろあると思いますが、rubyの側で、data02.memを読み込んで、Rプログラムの一部を書き出すという方法を採ります。すべてRで行うのも可能だと思いますが、私には少々敷居が高いので rubyを用います。
rubyの側で remake_factor() というメソッドを定義します。これは、次のように用います。
remake_factor("data02.mem", "dtf")
第１引数がファイル名、第２引数は Rプログラムにおけるデータフレームの名前です。
このメソッドの戻り値は下の２行です。
dtf$sex <- factor(dtf$sex, levels=c("m","f"), labels=c("男","女"))
{{br}}
dtf$view <- factor(dtf$view, levels=c("1","2","3"), labels=c("賛成","中立","反対"))
この２行をRプログラムの一部として挿入します。tjs05b.rb では erb を用いて挿入しています。
tjs05b.rb でやっていることは簡単なのですが、ちょっとごちゃごちゃしています。少しでも整理するため、rubyのメソッド定義を tjsf.rb に書き込んでおいて、それを取り込む形にしてあります。
erbを利用するサンプルにもなっているので、少し長いですが下に tjs05b.rb を全部掲げます。
−−−− tjs05b.rb ここから
#! ruby -Ks
# メモファイルに基づく数値・記号の置換＆erbの利用 (coding: Windows-31J)
require "erb"
require "rrxwin"
require "./tjsf"
##
rpro = <<EOS
dtf <- read.csv("<%= csv_file %>", header=T)
<%= remake_factor(memo_file, "dtf") %>
tbl <- table(性別=dtf$sex, 意見=dtf$view)
xx <- prop.table(addmargins(tbl,1), 1)*100
pct1 <- addmargins(xx, 2)
xx <- prop.table(addmargins(tbl,2), 2)*100
pct2 <- addmargins(xx, 1)
robj(list(tbl=addmargins(tbl), pct1=round(pct1,1), pct2=round(pct2,1)))
EOS
##
csv_file = "data02.csv"
memo_file = "data02.mem"
erb = ERB.new(rpro)
hs = Rrx.rexec(erb.result(binding))
xx = hs["robj1"]
keys = %w(tbl pct1 pct2)
titles = %w(度数分布 横の総数に対する構成比 縦の総数に対する構成比)
keys.each do |key|
puts titles.shift
print Rrx.ary2str(xx[key], "\t")
printf("\n")
end
−−−− tjs05b.rb ここまで
--------
!2. 欠損値のある単一回答の集計
回答者のID、sex(性別)、view(意見：賛成・中立・反対)が書かれたcsvファイル data03.csv を素材にします。
data02.csvと概ね同じですが、性別または意見の欄が空欄のケースがあります。いわば無回答があるわけです。
無回答がある場合の集計方法として、1)有効回答だけに着目する方法(無回答を除いた有効回答総数を100%とみなす)、2)無回答も１種の回答とみなす方法(無回答を含めた総数を100%とする)の２通りがあります。
どちらを採用するかは、調査項目の内容や、「無回答の背景にあると思われるもの」の判断の仕方によって違ってきます。
ということで、２通りの集計方法それぞれを取り上げてみます。
!!(1) 無回答を含めた集計
○ 要点１：欠損値を含む集計は「table(dtf$sex, dtf$view, exclude=NULL)」
○ 要点２：read.csv関数における欠損値の指定は na.strings オプション
有効回答のみに着目した集計、つまり無回答に全く言及しない集計は簡単です。tjs05.rbをそのまま使えます。
集計した度数分布は次のようになります。
|| ||賛成||中立||反対||Sum
||男||14||14||11||39
||女||9||17||20||46
||Sum||23||31||31||85
無回答が僅かで、考慮しなくてもいい旨を文章で説明していれば、上の表でも大丈夫かもしれません。
しかし、一般的には情報不足です。文章を読まずに表だけ拾い読みすることがあるので、やはり無回答がどれくらいかを表に盛り込む方が望ましいといえます。
Rプログラムでは、欠損値をNAという特別の値として扱うことが多いようです。NA は Not Available の意味。
table() 関数は、通常、NAを除いて集計しますが、「exclude=NULL」というオプションを指定すると、NAを含めた集計を行います。例えば次のようにします。
tbl <- table(性別=dtf$sex, 意見=dtf$view, exclude=NULL)
data03.csvの集計結果(度数分布)は下のとおりです。
|| ||賛成||中立||反対||NA||Sum
||男||14||14||11||4||43
||女||9||17||20||4||50
||NA||2||3||2||0||7
||Sum||25||34||33||8||100
「NA」の項をみると、無回答がどのくらいあるか把握できます。
以下にスクリプト tjs06.rb を示します。
−−−− tjs06.rb 抜粋ここから
rpro = <<EOS
dtf <- read.csv("data03.csv", header=T, na.strings="")
dtf$sex <- factor(dtf$sex, levels=c("m","f"), labels=c("男","女"))
dtf$view <- factor(dtf$view, levels=1:3, labels=c("賛成","中立","反対"))
tbl <- table(性別=dtf$sex, 意見=dtf$view, exclude=NULL)
xx <- prop.table(addmargins(tbl,1), 1)*100
pct1 <- addmargins(xx, 2)
xx <- prop.table(addmargins(tbl,2), 2)*100
pct2 <- addmargins(xx, 1)
robj(list(tbl=addmargins(tbl), pct1=round(pct1,1), pct2=round(pct2,1)))
EOS
##
hs = Rrx.rexec(rpro)
xx = hs["robj1"]
［後略］
−−−− tjs06.rb 抜粋ここまで
Rプログラムの read.csv() の最後の引数として「na.strings=""」を記述しています。これは、csvの空欄をNAにするとの指定です。
実は、この記述はなくてもいいのですが、念のため明示してみました。
もしcsvデータの中に NONE という４文字を書き込んでおいて、それを無効なものとして扱うのであれば、「na.strings="NONE"」とします。
性別欄は「m, f」の文字なので NONE にするとしても、意見の欄は「1, 2, 3」なので -999 をNAとして書き込みたい、ということであれば、次のようにします。
read.csv("data03.csv", header=T, na.strings=c("NONE",-999))
この場合、NONEがどこにあっても(性別欄以外のところであっても) NAになります。同じように、-999がどこにあっても(意見欄以外であっても) NAになるので注意して下さい。
!!［補足］ 行と列の項目名を変更する(「NA, Sum」→「無回答, 総数」)
tjs06.rb で作成した表には、NA, Sum という項目名がありました。これを「無回答」 「総数」に変更します。
xxが行列(tableやmatrix)である場合、colnames(xx) で列の項目名、rownames(xx) で行の項目名を得られます。これら両関数の戻り値は vector で、項目名が割り当てられていなければ NULL を返します。
列の項目名の中の NA を「無回答」に変更するには次のようにします。
colnames(xx)[is.na(colnames(xx))] <- "無回答"
また、列の項目名の中の「Sum」を「総数」に変更する場合は次のとおり。
colnames(xx)[colnames(xx)=="Sum"] <- "総数"
上の colnames を rownames にすれば、行の項目名を変更することになります。
これら変更の手続きをRプログラムの中に逐一 書き込むのは面倒なので、remake.item() という関数定義を tjsf.r に盛り込んで、それを呼び出すことにします。
xxが行列である場合、次のように用います。いくつかのパターンを示します。
xx <- remake.item(xx, c(NA,"Sum"), c("無回答","総数"))
xx <- remake.item(xx) # 「NA, Sum」→「無回答, 総数」
xx <- remake.item(xx, c("無効回答","総計")) # 「無効回答, 総計」に変更
基本的には３つの引数を指定します。第２引数は旧い名前群、第３引数は新しい名前群です。行と列の両方について、旧い名前群を新しい名前群に置き換えます。
第３引数を省略すると、第２引数が新しい名前群とみなされ、旧い名前群はデフォルト値の c(NA, "Sum") であるものとされます。
第２引数と第３引数の両方を省略すると、それぞれデフォルト値の c(NA, "Sum"), c("無回答", "総数") が指定されたものとみなされます。
--------
!!(2) 行列から行・列を取り出して引き算する(有効回答数の算出)
○ 要点：行列の5列目から4列目を引き算するには「xx[,5] - xx[,4]」
無回答を含む基本表(度数分布)を改めて示します。基本的に、tjs06.rbで作成される表と同じ内容です。
|| ||賛成||中立||反対||無回答||総数
||男||14||14||11||4||43
||女||9||17||20||4||50
||無回答||2||3||2||0||7
||総数||25||34||33||8||100
上の表をみても、有効回答がどれくらいかすぐには分かりません。
「総数」から「無回答」を引き算してやれば、有効回答を算出できるのでそれを試みます。とりあえず、縦方向に流れる列に着目します。
行列から特定の列を取り出すには xx[,4] とか xx[,5] とします。「無回答」が４列目、「総数」が５列目なので、その番号で取り出します。
また、項目名で取り出すこともできます。xx[,"無回答"], xx[,"総数"] と書くことができます。
つまり、有効回答は次のようにして求めることができます。２通りを示します。
yy <- xx[,5] - xx[,4]
yy <- xx[,"総数"] - xx[,"無回答"]
一方、縦方向の列ではなく、横方向の行を取り出して引き算する場合は次のとおりです。
行に着目すると、「無回答」が３番目、「総数」が４番目です。
yy <- xx[4,] - xx[3,]
yy <- xx["総数",] - xx["無回答",]
ここまで述べてきたノウハウを用いて、有効回答と無回答の両方を含む度数分布の表を作成してみます。
行にしろ列にしろ、「無回答」と「総数」が最後の２つの項目(つまり外側に位置する２つ)であることを利用して、行と列の取り出しを番号で行います。
また、行列の結合には cbind(), rbind() を使います。
以下、tjs07.rb です。
−−−− tjs07.rb 抜粋ここから
rpro = <<EOS
source("tjsf.r")
dtf <- read.csv("data03.csv", header=T, na.strings="")
dtf$sex <- factor(dtf$sex, levels=c("m","f"), labels=c("男","女"))
dtf$view <- factor(dtf$view, levels=1:3, labels=c("賛成","中立","反対"))
tbl <- table(性別=dtf$sex, 意見=dtf$view, exclude=NULL)
xx <- addmargins(tbl) # 行・列の両方に総数を付加
xx <- remake.item(xx) # 「NA, Sum」→「無回答, 総数」
n <- ncol(xx) # 列の個数
yy <- xx[,n] - xx[,n-1] # 有効回答数を算出
zz <- cbind(xx[,1:(n-2)], 有効回答=yy, xx[,(n-1):n])
xx <- zz
n <- nrow(xx) # 行の個数
yy <- xx[n,] - xx[n-1,] # 有効回答数を算出
zz <- rbind(xx[1:(n-2),], 有効回答=yy, xx[(n-1):n,])
robj(zz)
EOS
##
hs = Rrx.rexec(rpro)
zz = hs["robj1"]["self"]
print Rrx.ary2str(zz, "\t")
−−−− tjs07.rb 抜粋ここまで
cbind(), rbind() の第２引数を「有効回答=yy」としています。これを単に「yy」にすると、作成される表の項目名が「有効回答」にならず「yy」になります。
作成される表は次のとおりです。
|| ||賛成||中立||反対||有効回答||無回答||総数
||男||14||14||11||39||4||43
||女||9||17||20||46||4||50
||有効回答||23||31||31||85||8||93
||無回答||2||3||2||7||0||7
||総数||25||34||33||92||8||100
これで、有効回答がどれくらいか分かるようになりました。縦の有効回答と横の有効回答が交差する欄は 85 ですが、これは、性別と意見の両方とも分かる人の数です。
有効回答と無回答が交差する欄が何を意味するのか、ちょっと分かりにくい感じがありますが、性別が分かっているけれども意見が分からない人の数、その逆に、性別が分からないけれども意見が分かる人の数を示しています。
横の無回答と縦の無回答が交差する欄は、性別も意見も両方とも分からない人の数を示しますが、表ではそれが 0 になっています。
tjs07.rbのRプログラムの部分については、コメントを見ていただければ、どんなふうに処理を進めているのか分かると思います。
詳しい仕様は省略しますが、cbind() は、行列に列を追加する関数です。a, b が行列やベクトルである場合、「x<-cbind(a,b)」とすれば、aにbを追加した結果がxに入ります。xは、列を追加した結果なので、aよりも横に拡がったものになります。cbind(a,b,c) のように３つ以上の引数を与えることもできます。
rbind() は、行を追加するものです。追加した結果は、当然ながら縦に拡がります。
--------
!!(3) 総数に対する比率(パーセント)の表を作成する
○ 要点：１行目について最右端の列で割り算した結果(パーセント)を得るには{{br}}
「xx[1,] / xx[1,ncol(xx)] * 100」
tjs07.rbで作成した表は、度数分布の表です。これを基にして、各欄の値が「総数」に占める比率(パーセント)を求めてみます。そうすれば、有効回答率等が分かります。
ちなみに、有効回答を100%とする比率(構成比)は、tjs05.rbが出力する表で把握できます。
構成比を求めるのに、前のサンプルでは prop.table() を使いました。tjs03.rb がそれです。
しかし、今回は prop.table を利用しません。「有効回答」が一種の小計になっているため、「総数」が各欄を足し算したものにならないからです。
例えば、「男」に着目して、「賛成・中立・反対・有効回答・無回答」を素直に足し算すると、82 になります。これは「総数」の 43 よりずっと大きい数です。「有効回答」の分だけ多くなっています。
そこで、各欄を「総数」で割り算して比率を求めることにします。
とりあえず第１行目の「男」に着目すると、「総数」が最右端の６列目なので、「総数」を含めた各欄の比率(パーセント)は次のようにして算出できます。
male <- xx[1,] / xx[1,6] * 100
こうすると、「賛成・中立・反対・有効回答・無回答・総数」の各パーセントの値が male に入ります。
「男」だけでなく「女・有効回答・無回答・総数」の５行分について比率を求めようとすれば、次のようになります。基本表は xx に入っているものとします。
nc <- ncol(xx) # 最右端の列の番号
for (i in 1:nrow(xx)) {
if (xx[i,nc] == 0) next # 0による割り算は避ける
xx[i,] <- xx[i,] / xx[i,nc] * 100 # i行目を比率に入れ換え
}
上の処理を経ると、度数が記録されていたxxの中身が、比率(パーセント)に入れ換えられます。
tjs07.rbで得られる基本表に対し、上の処理を施して比率を得るようにしたのが tjs08.rb です。
参考まで掲げておきます。
−−−− tjs08.rb 抜粋ここから
rpro = <<EOS
source("tjsf.r")
dtf <- read.csv("data03.csv", header=T, na.strings="")
dtf$sex <- factor(dtf$sex, levels=c("m","f"), labels=c("男","女"))
dtf$view <- factor(dtf$view, levels=1:3, labels=c("賛成","中立","反対"))
tbl <- table(性別=dtf$sex, 意見=dtf$view, exclude=NULL)
xx <- addmargins(tbl) # 行・列の両方に総数を付加
xx <- remake.item(xx) # 「NA, Sum」→「無回答, 総数」
n <- ncol(xx) # 列の個数
yy <- xx[,n] - xx[,n-1] # 有効回答数を算出
zz <- cbind(xx[,1:(n-2)], 有効回答=yy, xx[,(n-1):n])
xx <- zz
n <- nrow(xx) # 行の個数
yy <- xx[n,] - xx[n-1,] # 有効回答数を算出
tbl <- rbind(xx[1:(n-2),], 有効回答=yy, xx[(n-1):n,])
nc <- ncol(tbl)
nr <- nrow(tbl)
pct1 <- tbl
for (i in 1:nr) {
if (pct1[i,nc] == 0) next
pct1[i,] <- pct1[i,] / pct1[i,nc] * 100
}
pct2 <- tbl
for (i in 1:nc) {
if (pct2[nr,i] == 0) next
pct2[,i] <- pct2[,i] / pct2[nr,i] * 100
}
robj(list(tbl=tbl, pct1=round(pct1,1), pct2=round(pct2,1)))
EOS
##
hs = Rrx.rexec(rpro)
xx = hs["robj1"]
［後略］
−−−− tjs08.rb 抜粋ここまで
--------
以上で「単純集計に関する覚え書き／単一回答の集計」をおわりにします。
Copyright (C) T. Yoshiizumi, 2013 All rights reserved.