社会調査データを扱う場合に実際に面倒なのは,データ分析や多変量解析のくだりではなく,むしろその準備作業としてのデータの整備や加工の部分である。収集した情報(調査票の回答選択肢番号)から別の数値に割り当てるリコードや,情報が得られなかった場合の欠損値の処理などがきちんと出来ないと,そもそも適切な集計や分析が行えない。しかしこれらの作業は取り分け統計ソフトによって異なる部分であり,特にRの場合にはSPSSなどとやや発想を変えなければうまく使えない事も多い。
心理統計学では変数の種類(尺度水準)を,名義尺度,順序尺度,間隔尺度,比例尺度に分ける事が多い。名義尺度や順序尺度の情報を数値で表現する事は多いが,足し算や引き算などの演算が出来るのは間隔尺度と比例尺度である。
他方,SPSSなどの統計ソフトをはじめ,コンピュータソフトは一般的に,数値型/文字型の様に情報を区別し,数値型と認識されたものには加減乗除などの演算を行う。数字でさえあれば,それが名義尺度であるか順序尺度であるか間隔尺度以上であるかは問題としない。
コンピュータソフトを用いて統計分析を行う場合,変数の種類或いは尺度水準を正しく理解して区別しないと,気付かないうちに誤った分析をしてしまう事になるので,この点についてよく注意しなければならない。
Rでは取り敢えず,数値型numeric,文字型character,要因型factorの3つの区別に注意する事にしよう。
例えば,社会調査で性別を訪ねて,男性=1,女性=2,その他=3として回答を収集したとしよう。これは,文字で表現しようが数字で表現しようが,名義尺度である。
10人の回答を数値でヴェクトルq01に格納したとする。
q01 <- c(1, 2, 3, 2, 2, 1, 1, 2, 1, 3)
Rの関数 mean( ) や summary( ) をq01に適用すると,算術平均や中央値などを出力するが,これらはいずれも名義尺度について求める事の出来ない統計量である。 このヴェクトルの型を調べる為に class( ) と云う関数を適用すると “numeric” と表示される。numeric,数値型変数と認識されている限り,平均でも分散でも計算されてしまう。
q01
[1] 1 2 3 2 2 1 1 2 1 3
mean(q01)
[1] 1.8
summary(q01)
Min. 1st Qu. Median Mean 3rd Qu. Max.
1.0 1.0 2.0 1.8 2.0 3.0
class(q01)
[1] "numeric"
今度は,同じ回答が,文字のかたちでヴェクトルq02に格納されているとしよう。
q02 <- c("男性", "女性", "その他", "女性", "女性", "男性", "男性", "女性", "男性", "その他")
こうすると mean( ) は計算が拒否され,summary( ) もデータの個数(Length)以外は分からない。
class( ) とすると,“character”,文字列である事が分かる(summary( )の結果にも表示されている)。
q02
[1] "男性" "女性" "その他" "女性" "女性" "男性" "男性"
[8] "女性" "男性" "その他"
mean(q02)
[1] NA
summary(q02)
Length Class Mode
10 character character
class(q02)
[1] "character"
Rにはこれ以外に要因型factorと云う変数があり,データ分析では数値型と並んでこれを使う事になる。
単なる文字型ヴェクトルのq02から,要因型のヴェクトルq03を作成してみよう。
q03 <- factor(q02)
summary(q03)
その他 女性 男性
2 4 4
class(q03)
[1] "factor"
算術平均は勿論計算されないが,summary( ) を使うと度数分布表が出力される点でcharacterと異なる。class(q03) とすると“factor”と出力される。q02; q03としてただ単に中身を表示させるだけでも,factorの方だけ「Levels: その他 女性 男性」と云う行が出力されている。
q04 <- factor(q02, levels=c("男性", "女性", "その他"))
summary(q04)
男性 女性 その他
4 4 2
要因型変数は名義尺度として扱われる。名義尺度や順序尺度のカテゴリカル変数を分析する際には,きちんと要因型の“factor”になっているかどうかを確認し,“numeric”になっていれば上の様に“factorに”変換した新変数を作成してから分析に用いよう。
数値型変数から要因型変数を作成する事も出来る。
q05 <- factor(q01, levels = 1:3, labels = c("1 男性", "2 女性", "3 その他"))
summary(q05); class(q05)
1 男性 2 女性 3 その他
4 4 2
[1] "factor"
社会調査では,情報を得ようとしても得られない場合が多々ある。得られなかった情報を「欠損」や「欠測」,missing value, nonresponse(item nonresponse, unit nonresponse) と呼ぶ。
MS-Excelやcsvでデータ・ファイルを作成する際,欠損値には,予め決めておいたルールに従って9や99を入力したり,単に空欄にしておいたりする。9や99などの数値が入力されている場合は,それを読み込んで作成されたRのデータ・フレイムでも9や99として他の数値と同様に扱われる。特にそれが欠損値を意味すると云う事は考慮されないので注意が必要である。csvファイルで空欄になっているセルについては,Rでの欠損値を意味する NA としてデータ・フレイムに読み込まれる。
9や99として読み込まれると,うっかり有効なデータだとして計算に含まれたりするので注意が必要であり,NAとして読み込まれると,うっかり有効扱いにされる事は無いが,逆に関数が意図した計算結果を返さない事があったりして面倒である。
次の様に,最後の欠損値の部分だけが異なる二つのヴェクトルで試してみよう。summary(q06a)とsummary(q06b)の結果を見比べると,q06aは其の儘で使用出来ない事が分かる。
q06a <- c(12, 12, 14, 18, 16, 16, 9, 16, 12, 16, 99)
q06b <- c(12, 12, 14, 18, 16, 16, 9, 16, 12, 16, NA)
summary(q06a); summary(q06b)
Min. 1st Qu. Median Mean 3rd Qu. Max.
9.00 12.00 16.00 21.82 16.00 99.00
Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
9.0 12.0 15.0 14.1 16.0 18.0 1
しかし,q06bは,summary( )関数では問題を生じないが,mean( ) や sd( ) と云った基本的な関数で,其の儘では計算を拒否されてしまう。いずれも,結果が NA で返ってきてしまう。
これを防ぐには,mean(q06b, na.rm=T) や sd(q06b, na.rm=T) の様に,na.rm=T と云うオプションを追加しなければならない。
mean(q06b); sd(q06b); var(q06b)
[1] NA
[1] NA
[1] NA
mean(q06b, na.rm=T); sd(q06b, na.rm=T); var(q06b, na.rm=T)
[1] 14.1
[1] 2.766867
[1] 7.655556
また,度数分布表を作成しようと table(q06b) とすると今度はNAが自動的に除外されて出力されるが,これでは NA が存在するのかどうか分からない。データやファイルのチェック作業などでは,欠損値を含んで度数分布表や分割表(クロス表)を作成しなければならないので不便である。この場合は table(q06b, useNA=“always”) 又は table(q06b, useNA=“ifany”) としなければならない。useNA=のオプションを省略した場合は自動的に useNA=“no” として扱われているのである。
table(q06b)
q06b
9 12 14 16 18
1 3 1 4 1
table(q06b, useNA="always")
q06b
9 12 14 16 18 <NA>
1 3 1 4 1 1
table(q06b, useNA="ifany")
q06b
9 12 14 16 18 <NA>
1 3 1 4 1 1
table(q06b, useNA="no")
q06b
9 12 14 16 18
1 3 1 4 1
その他,相関係数を求める関数 cor( ) も,NAが含まれていると計算してくれないので,cor(x, y , use=“complete”) の様に,use=“complete”もしくはuse=“pairwise”のオプションを指定しなければならない。
同じ相関係数でも,検定や区間推定を同時に行う cor.test( ) は自動的にNAを除外して計算する。
q07b <- c(400, NA, 600, 1000, 1500, 700, 200, 500, 800, 900, NA)
cor(q06b, q07b)
[1] NA
cor(q06b, q07b , use="complete")
[1] 0.6499337
cor.test(q06b, q07b)
Pearson's product-moment correlation
data: q06b and q07b
t = 2.2626, df = 7, p-value = 0.05811
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
-0.02496288 0.91787036
sample estimates:
cor
0.6499337
データ集計や分析で何故か計算結果が NA として返される場合には,データの中にNAが含まれているのに適切にそれを除外していない可能性を考えて見よう。また table( ) を使用する時には,useNA=“always”のオプションを付けた場合と見比べて,NAの有無を確認する様にしよう。
欠損値を処理するオプションの指定の仕方は関数によって異なっている事もあるので,知っているオプションでうまく行かない場合にはウェブなどで調べてみよう。
社会調査では回答選択肢の数字を其の儘データファイルに記録するのが通例だが,回答選択肢番号が分析可能な数量変数になるとは限らない。
例えば,最終学歴schoolが「中学=1,高校=2,短期大学=3,四年制大学=4」として記録されている場合,これを教育年数(それぞれ,9,12,14,16)として分析に用いたい事が多い。
仮に10人分のデータ school <- c(3, 3, 2, 2, 4, 2, 4, 3, 1, 4, 4) であるとしよう。
SPSSであれば recode と云うコマンドを使うが,Rに同じものは無い。どんなコンピュータ言語にもあるif文と云うのを使う方法もあるが,schoolはあくまで一つの数字ではなくヴェクトルである事をきちんと理解していないと,意図した動作をしない。
「schoolが1ならばeduを9にする」を其の儘コマンドに表現すると, if (school == 1) edu <- 9 となりそうだが,school == 1 はschoolが一つの数字である事を想定している表現となっておりうまく行かない。表示されるエラーメッセイジと edu の内容を確認すれば分かる。
school <- c(3, 3, 2, 2, 4, 2, 4, 3, 1, 4, 4)
if (school == 1) edu <- 9
edu
Error in eval(expr, envir, enclos): オブジェクト 'edu' がありません
school自体はヴェクトルであり,その中の要素が1に等しいか否かを評価したいので,以下の様にすると意図した結果になる。
edu <- NULL
edu[school == 1] <- 9
edu
[1] NA NA NA NA NA NA NA NA 9 NA NA
[ ] はヴェクトル中の要素を指定する為のものであった。この [ ] の中に書いてあるものは要素についての条件であると解釈され,“schoolの値が1であるケース(行)についてはeduの値に9を代入せよ”と云う意味になる。但しこの場合,既にeduと云うヴェクトルが存在している事を前提とした表現になるので,先にeduを空箱(NULL)として作成しておく。作成されたeduを確認すると一ヵ所だけ9になっており,他はNAが代入されている。この方式で意図した変換は可能だが,schoolの全ての値について式を書かなければならない。
もっと効率的な方法として,次の様な方法がある。
c(9, 12, 14, 16)はヴェクトルである。c(9, 12, 14, 16)の要素を指定したい場合(例えば2番目の要素である12),c(9, 12, 14, 16)[2]などとする事になる。
要素は複数を指定する事が出来るので,c(9, 12, 14, 16)[c(2, 3)]などとする事が出来る([2, 3]と指定すると二次元ヴェクトルの2行3列目を指示している事になるのでエラーになる)。
schoolも一次元ヴェクトルであるので,こうしたヴェクトルの要素指定を応用してSPSSのrecode同様に(或いはより簡単に)新変数を作成する事が出来るのである。
c(9, 12, 14, 16)
[1] 9 12 14 16
c(9, 12, 14, 16)[2]
[1] 12
edu <- c(9, 12, 14, 16)[school]
edu
[1] 14 14 12 12 16 12 16 14 9 16 16
新変数を作成したら,必ず以下の様に正しく作成されている事を確認する。
table(school, edu, useNA = "always")
edu
school 9 12 14 16 <NA>
1 1 0 0 0 0
2 0 3 0 0 0
3 0 0 3 0 0
4 0 0 0 4 0
<NA> 0 0 0 0 0
ここで,最初のデータに欠損値の意味で99が記録されている場合を考えよう。
school <- c(3, 3, 2, 2, 4, 2, 4, 3, 1, 4, 99)
c(9, 12, 14, 16)にも99番目の要素を追加するのは(簡単に出来る事は出来るが)面倒である。また,上の欠損値の箇所で述べた様に,99はNAとして読み込んだ方が良い。
実は先の edu <- c(9, 12, 14, 16)[school] を其の儘用いると,99をNAにする変換も一緒に出来る。ヴェクトルc(9, 12, 14, 16)の99番目の要素を指定しようとすると,そんな要素は存在しないのでNAになるのである。よって,欠損値に変換したい値については,わざわざ準備する必要が無い。
edu <- c(9, 12, 14, 16)[school]
edu
[1] 14 14 12 12 16 12 16 14 9 16 NA
table(school, edu, useNA = "always")
edu
school 9 12 14 16 <NA>
1 1 0 0 0 0
2 0 3 0 0 0
3 0 0 3 0 0
4 0 0 0 3 0
99 0 0 0 0 1
<NA> 0 0 0 0 0
本人学歴school,配偶者学歴sschool,父学歴fschool,母学歴mschoolに全て同じ変換を施したい場合,変換に用いるヴェクトルに名前を付けて使用するとすっきりして見易い。
x <- c(9, 12, 14, 16)
edu <- x[school]; sedu <- x[sschool]; fedu <- x[fschool]; medu <- x[mschool]
年収などのデータも,本人年収,配偶者年収,世帯年収などの情報を取得していたりするので,こうした変換が便利である。以下は,データ・フレイムdata01中の本人年収q24,配偶者年収q25,世帯年収q26の選択肢番号情報を,年収額を表すincome,sincome,fincomeに変換する例である。選択肢番号1から16には年収カテゴリが割り当てられており,17は「わからない」,18は「答えたくない」となっているものとする。17,18は欠損値にしたいので,1~16に対応する金額を指定するだけで良い。
x_inc <- c(0, 25, 75, 125, 200, 300, 400, 500, 600, 700, 800, 925, 1125, 1375, 1750, 2250)
data01$income <- x_inc[data01$q24]
data01$sincome <- x_inc[data01$q25]
data01$fincome <- x_inc[data01$q26]
性別はよく「男性=1,女性=2」(調査によっては更に「その他=3」等がある)と記録されているが,俗に「男性ダミー」(或いは女性ダミー)と呼ばれる変数に変換したい事がしばしばある。男性ダミー変数とは,男性であれば1,そうでなければ0の値をとる変数の事である。変数sexは1(男性),2(女性),3(その他)の値を取るとすると,男性ダミーmaleは次の様に作成出来る。
male <- c(1, 0, 0)[sex]
女性ダミーなら female <- c(0, 1, 0)[sex] である。
上のshool変数を,初等学歴(義務教育)/中等学歴(高校)/高等学歴(短大以上)の3区分の名義尺度gakuに変換するには次の様にすれば良い。
gaku0 <- c(1, 2, 3, 3)[school]
gaku <- factor(gaku0, levels=1:3, labels=c("1 初等", "2 中等", "3 高等"))
gaku
[1] 3 高等 3 高等 2 中等 2 中等 3 高等 2 中等 3 高等 3 高等 1 初等 3 高等
[11] <NA>
Levels: 1 初等 2 中等 3 高等
summary(gaku)
1 初等 2 中等 3 高等 NA's
1 3 6 1