進化学実習 2023 牧野研 東北大学
東北大学 生命科学研究科 進化ゲノミクス分野 牧野研 特任助教
- 導入: データ解析の全体像。Rの基本。
- データの可視化、レポート作成。
- データ構造の処理1: 抽出、集約など。
- データ構造の処理2: 結合、変形など。
- データ内容の処理: 数値、文字列など。
- データ入力、データ解釈
- 統計モデリング1: 確率分布、尤度
- 統計モデリング2: 一般化線形モデル
- 発表会
2023-04-12 東北大学 理学部生物学科 進化学実習
https://heavywatal.github.io/slides/tohoku2023r/
https://heavywatal.github.io/slides/tohoku2023r/
データ解析のおおまかな流れ
- コンピュータ環境の整備 ✅
- データの取得、読み込み 👈 今回はここから
- 探索的データ解析
- 前処理、加工 ✅ 完全に理解した
- 可視化、仮説生成 ✅ 完全に理解した
- 統計解析、仮説検証 ⬜ 今日の終盤と次回
- 報告、発表 ✅ Quarto楽しい
data.frameの読み書き
-
readxlパッケージを使えば
.xlsxファイルも読める、けど -
カンマ区切り(CSV)とかタブ区切り(TSV)のテキストが無難。
-
R標準の関数ではなく readr パッケージを使う。
-
ファイル名は作業ディレクトリからの相対パスで指定。
readr::write_tsv(iris, "data/iris.tsv") # 書き出し iris2 = readr::read_tsv("data/iris.tsv") # 読み込み
あれれー、エラーが出る？
Error: Cannot open file for writing:
* 'data/iris.tsv'
ファイル読み書きでありがちなエラー
- 読みたいファイル名が間違っている。
- 読みたいファイルが作業ディレクトリとは違うところにある。
- 書き込み先のディレクトリが存在しない。→ 予め作っておこう
冷静に、現在の作業ディレクトリとその中身を確認しよう:
getwd() # GET Working Directory
fs::dir_ls(".") # List files in "." (here)
fs::dir_ls("data") # List files in "./data"
fs::dir_create("data") # Create directory
よくあるエラー集 (石川由希さん@名古屋大) を読んでおきましょう。
🔰 これまでに登場したデータを data/ 以下に書き出してみよう。
🔰 それらを別の名前で読み込み直してみよう。
🔰 最終課題予告
読み込み → 前処理 → 可視化 → 余裕でクリア、、、となるか？
実演: e-Stat 国勢調査 CSV
2020年 → 小地域 → 年齢（５歳階級、４区分）別、男女別人口 → 宮城県
うまくすれば次のような人口ピラミッドを描けるはずのファイルだが…
オープンデータを入手できても簡単には読み込めない
初手エラー:
infile = "tblT001082C04.txt"
readr::read_csv(infile)
Error in nchar(x, keepNA = FALSE): invalid multibyte string, element 2
RStudioでテキストとして開いてみると文字化け:
KEY_CODE,HYOSYO,CITYNAME,NAME,HTKSYORI,HTKSAKI,GASSAN,T001082001,T001082002,T001082003,T001082004,T001082005,T001082006,T001082007,T001082008,T001082009,T001082010,T001082011,T001082012,T001082013,T001082014,T001082015,T001082016,T001082017,T001082018,T001082019,T001082020,T001082021,T001082022,T001082023,T001082024,T001082025,T001082026,T001082027,T001082028,T001082029,T001082030,T001082031,T001082032,T001082033,T001082034,T001082035,T001082036,T001082037,T001082038,T001082039,T001082040,T001082041,T001082042,T001082043,T001082044,T001082045,T001082046,T001082047,T001082048,T001082049,T001082050,T001082051,T001082052,T001082053,T001082054,T001082055,T001082056,T001082057,T001082058,T001082059,T001082060
,,,,,,,<91><8d><90><94><81>A<94>N<97><ee><81>u<95>s<8f>ځv<8a>܂<de>,<91><8d><90><94><82>O<81>`<82>S<8d><ce>,<91><8d><90><94><82>T<81>`<82>X<8d><ce>,<91><8d><90><94><82>P<82>O<81>`<82>P<82>S<8d><ce>,<91><8d><90><94><82>P<82>T<81>`<82>P<82>X<8d><ce>,<91><8d><90><94><82>Q<82>O<81>`<82>Q<82>S<8d><ce>,<91><8d><90><94><82>Q<82>T<81>`<82>Q<82>X<8d><ce>,<91><8d><90><94><82>R<82>O<81>`<82>R<82>S<8d><ce>,<91><8d><90><94><82>R<82>T<81>`<82>R<82>X<8d><ce>,<91><8d><90><94><82>S<82>O<81>`<82>S<82>S<8d><ce>,<91><8d><90><94><82>S<82>T<81>`<82>S<82>X<8d><ce>,<91><8d><90><94><82>T<82>O<81>`<82>T<82>S<8d><ce>,<91><8d><90><94><82>T<82>T<81>`<82>T<82>X<8d><ce>,<91><8d><90><94><82>U<82>O<81>`<82>U<82>S<8d><ce>,<91><8d><90><94><82>U<82>T<81>`<82>U<82>X<8d><ce>,<91><8d><90><94><82>V<82>O<81>`<82>V<82>S<8d><ce>,<91><8d><90><94><82>P<82>T<8d>Ζ<a2><96><9e>,<91><8d><90><94><82>P<82>T<81>`<82>U<82>S<8d><ce>,<91><8d><90><94><82>U<82>T<8d>Έȏ<e3>,<91><8d><90><94><82>V<82>T<8d>Έȏ<e3>,<92>j<82>̑<8d><90><94><81>A<94>N<97><ee><81>u<95>s<8f>ځv<8a>܂<de>,<92>j<82>O<81>`<82>S<8d><ce>,<92>j<82>T<81>`<82>X<8d><ce>,<92>j<82>P<82>O<81>`<82>P<82>S<8d><ce>,<92>j<82>P<82>T<81>`<82>P<82>X<8d><ce>,<92>j<82>Q<82>O<81>`<82>Q<82>S<8d><ce>,<92>j<82>Q<82>T<81>`<82>Q<82>X<8d><ce>,<92>j<82>R<82>O<81>`<82>R<82>S<8d><ce>,<92>j<82>R<82>T<81>`<82>R<82>X<8d><ce>,<92>j<82>S<82>O<81>`<82>S<82>S<8d><ce>,<92>j<82>S<82>T<81>`<82>S<82>X<8d><ce>,<92>j<82>T<82>O<81>`<82>T<82>S<8d><ce>,<92>j<82>T<82>T<81>`<82>T<82>X<8d><ce>,<92>j<82>U<82>O<81>`<82>U<82>S<8d><ce>,<92>j<82>U<82>T<81>`<82>U<82>X<8d><ce>,<92>j<82>V<82>O<81>`<82>V<82>S<8d><ce>,<92>j<82>P<82>T<8d>Ζ<a2><96><9e>,<92>j<82>P<82>T<81>`<82>U<82>S<8d><ce>,<92>j<82>U<82>T<8d>Έȏ<e3>,<92>j<82>V<82>T<8d>Έȏ<e3>,<8f><97><82>̑<8d><90><94><81>A<94>N<97><ee><81>u<95>s<8f>ځv<8a>܂<de>,<8f><97><82>O<81>`<82>S<8d><ce>,<8f><97><82>T<81>`<82>X<8d><ce>,<8f><97><82>P<82>O<81>`<82>P<82>S<8d><ce>,<8f><97><82>P<82>T<81>`<82>P<82>X<8d><ce>,<8f><97><82>Q<82>O<81>`<82>Q<82>S<8d><ce>,<8f><97><82>Q<82>T<81>`<82>Q<82>X<8d><ce>,<8f><97><82>R<82>O<81>`<82>R<82>S<8d><ce>,<8f><97><82>R<82>T<81>`<82>R<82>X<8d><ce>,<8f><97><82>S<82>O<81>`<82>S<82>S<8d><ce>,<8f><97><82>S<82>T<81>`<82>S<82>X<8d><ce>,<8f><97><82>T<82>O<81>`<82>T<82>S<8d><ce>,<8f><97><82>T<82>T<81>`<82>T<82>X<8d><ce>,<8f><97><82>U<82>O<81>`<82>U<82>S<8d><ce>,<8f><97><82>U<82>T<81>`<82>U<82>X<8d><ce>,<8f><97><82>V<82>O<81>`<82>V<82>S<8d><ce>,<8f><97><82>P<82>T<8d>Ζ<a2><96><9e>,<8f><97><82>P<82>T<81>`<82>U<82>S<8d><ce>,<8f><97><82>U<82>T<8d>Έȏ<e3>,<8f><97><82>V<82>T<8d>Έȏ<e3>
04101,1,<90><e5><91><e4><8e>s<90>t<8b><e6>,,0,,,311590,10231,11633,11838,15944,23772,17838,17677,19028,21113,23299,20544,18295,16702,17137,17700,33702,194212,69969,35132,150932,5350,6087,6066,8321,12579,8669,8600,9241,10232,11412,10317,9031,7893,8010,8082,17503,96295,29477,13385,160658,4881,5546,5772,7623,11193,9169,9077,9787,10881,11887,10227,9264,8809,9127,9618,16199,97917,40492,21747
041010010,2,<90><e5><91><e4><8e>s<90>t<8b><e6>,<90>t<92><ac>,0,,,649,16,15,17,23,53,62,49,40,40,40,45,33,28,38,38,48,413,143,67,307,8,10,7,10,26,30,27,21,18,14,25,16,11,15,16,25,198,60,29,342,8,5,10,13,27,32,22,19,22,26,20,17,17,23,22,23,215,83,38
041010020,2,<90><e5><91><e4><8e>s<90>t<8b><e6>,<82><a0><82><af><82>ڂ̒<ac>,0,,,741,23,18,13,26,32,55,48,42,60,51,48,47,38,55,43,54,447,209,111,365,12,10,4,15,16,22,28,19,39,27,25,16,23,24,21,26,230,89,44,376,11,8,9,11,16,33,20,23,21,24,23,31,15,31,22,28,217,120,67
メニューから “File → Reopen with Encoding…” を選択。
たいてい、古い日本語ファイルにありがちなSHIFT-JISが正解。
(もっと古いファイルならEUC-JPとかの場合もある。)
文字コードを指定して読み込み
よく見ると、8列目以降の列名が2行にまたがっている:
sjis = readr::locale(encoding = "SHIFT-JIS")
readr::read_csv(infile, locale = sjis)
KEY_CODE HYOSYO CITYNAME NAME HTKSYORI HTKSAKI GASSAN T001082001 T001082002 T001082003 T001082004 T001082005 T001082006 T001082007 T001082008 T001082009 T001082010 T001082011 T001082012 T001082013 T001082014 T001082015 T001082016 T001082017 T001082018 T001082019 T001082020 T001082021 T001082022 T001082023 T001082024 T001082025 T001082026 T001082027 T001082028 T001082029 T001082030 T001082031 T001082032 T001082033 T001082034 T001082035 T001082036 T001082037 T001082038 T001082039 T001082040 T001082041 T001082042 T001082043 T001082044 T001082045 T001082046 T001082047 T001082048 T001082049 T001082050 T001082051 T001082052 T001082053 T001082054 T001082055 T001082056 T001082057 T001082058 T001082059 T001082060
1 <NA> NA <NA> <NA> NA <NA> <NA> 総数、年齢「不詳」含む 総数０〜４歳 総数５〜９歳 総数１０〜１４歳 総数１５〜１９歳 総数２０〜２４歳 総数２５〜２９歳 総数３０〜３４歳 総数３５〜３９歳 総数４０〜４４歳 総数４５〜４９歳 総数５０〜５４歳 総数５５〜５９歳 総数６０〜６４歳 総数６５〜６９歳 総数７０〜７４歳 総数１５歳未満 総数１５〜６４歳 総数６５歳以上 総数７５歳以上 男の総数、年齢「不詳」含む 男０〜４歳 男５〜９歳 男１０〜１４歳 男１５〜１９歳 男２０〜２４歳 男２５〜２９歳 男３０〜３４歳 男３５〜３９歳 男４０〜４４歳 男４５〜４９歳 男５０〜５４歳 男５５〜５９歳 男６０〜６４歳 男６５〜６９歳 男７０〜７４歳 男１５歳未満 男１５〜６４歳 男６５歳以上 男７５歳以上 女の総数、年齢「不詳」含む 女０〜４歳 女５〜９歳 女１０〜１４歳 女１５〜１９歳 女２０〜２４歳 女２５〜２９歳 女３０〜３４歳 女３５〜３９歳 女４０〜４４歳 女４５〜４９歳 女５０〜５４歳 女５５〜５９歳 女６０〜６４歳 女６５〜６９歳 女７０〜７４歳 女１５歳未満 女１５〜６４歳 女６５歳以上 女７５歳以上
2 04101 1 仙台市青葉区 <NA> 0 <NA> <NA> 311590 10231 11633 11838 15944 23772 17838 17677 19028 21113 23299 20544 18295 16702 17137 17700 33702 194212 69969 35132 150932 5350 6087 6066 8321 12579 8669 8600 9241 10232 11412 10317 9031 7893 8010 8082 17503 96295 29477 13385 160658 4881 5546 5772 7623 11193 9169 9077 9787 10881 11887 10227 9264 8809 9127 9618 16199 97917 40492 21747
3 041010010 2 仙台市青葉区 青葉町 0 <NA> <NA> 649 16 15 17 23 53 62 49 40 40 40 45 33 28 38 38 48 413 143 67 307 8 10 7 10 26 30 27 21 18 14 25 16 11 15 16 25 198 60 29 342 8 5 10 13 27 32 22 19 22 26 20 17 17 23 22 23 215 83 38
4 041010020 2 仙台市青葉区 あけぼの町 0 <NA> <NA> 741 23 18 13 26 32 55 48 42 60 51 48 47 38 55 43 54 447 209 111 365 12 10 4 15 16 22 28 19 39 27 25 16 23 24 21 26 230 89 44 376 11 8 9 11 16 33 20 23 21 24 23 31 15 31 22 28 217 120 67
--
5941 04606004015 4 南三陸町 歌津字石浜 0 <NA> <NA> 295 9 7 11 6 6 12 6 18 16 19 21 36 18 30 26 27 158 110 54 146 5 5 5 3 5 8 2 6 8 10 14 16 8 11 20 15 80 51 20 149 4 2 6 3 1 4 4 12 8 9 7 20 10 19 6 12 78 59 34
5942 04606004016 4 南三陸町 歌津字田の浦 0 <NA> <NA> 144 5 2 5 3 5 7 5 7 6 5 15 12 17 14 4 12 82 50 32 66 - 1 3 1 3 4 3 4 3 2 6 5 9 9 4 4 40 22 9 78 5 1 2 2 2 3 2 3 3 3 9 7 8 5 - 8 42 28 23
5943 04606004017 4 南三陸町 歌津字草木沢 0 <NA> <NA> 457 21 16 18 16 12 13 20 18 18 41 34 28 37 46 43 55 237 165 76 234 11 6 8 10 11 5 9 10 11 22 18 15 17 19 29 25 128 81 33 223 10 10 10 6 1 8 11 8 7 19 16 13 20 27 14 30 109 84 43
5944 04606004018 4 南三陸町 歌津字伊里前 0 <NA> <NA> - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
左側と右側を分けて読み、それぞれ削ってからbind
なんか数値の列のはずなのに - という文字が混じっている。(中には X も):
dfL = readr::read_csv(infile, locale = sjis, col_select = seq(1, 7)) |>
dplyr::slice(-1)
dfR = readr::read_csv(infile, locale = sjis, col_select = -seq(1, 7),
skip = 1L)
raw_miyagi = dplyr::bind_cols(dfL, dfR) |> print()
KEY_CODE HYOSYO CITYNAME NAME HTKSYORI HTKSAKI GASSAN 総数、年齢「不詳」含む 総数０〜４歳 総数５〜９歳 総数１０〜１４歳 総数１５〜１９歳 総数２０〜２４歳 総数２５〜２９歳 総数３０〜３４歳 総数３５〜３９歳 総数４０〜４４歳 総数４５〜４９歳 総数５０〜５４歳 総数５５〜５９歳 総数６０〜６４歳 総数６５〜６９歳 総数７０〜７４歳 総数１５歳未満 総数１５〜６４歳 総数６５歳以上 総数７５歳以上 男の総数、年齢「不詳」含む 男０〜４歳 男５〜９歳 男１０〜１４歳 男１５〜１９歳 男２０〜２４歳 男２５〜２９歳 男３０〜３４歳 男３５〜３９歳 男４０〜４４歳 男４５〜４９歳 男５０〜５４歳 男５５〜５９歳 男６０〜６４歳 男６５〜６９歳 男７０〜７４歳 男１５歳未満 男１５〜６４歳 男６５歳以上 男７５歳以上 女の総数、年齢「不詳」含む 女０〜４歳 女５〜９歳 女１０〜１４歳 女１５〜１９歳 女２０〜２４歳 女２５〜２９歳 女３０〜３４歳 女３５〜３９歳 女４０〜４４歳 女４５〜４９歳 女５０〜５４歳 女５５〜５９歳 女６０〜６４歳 女６５〜６９歳 女７０〜７４歳 女１５歳未満 女１５〜６４歳 女６５歳以上 女７５歳以上
1 04101 1 仙台市青葉区 <NA> 0 <NA> <NA> 311590 10231 11633 11838 15944 23772 17838 17677 19028 21113 23299 20544 18295 16702 17137 17700 33702 194212 69969 35132 150932 5350 6087 6066 8321 12579 8669 8600 9241 10232 11412 10317 9031 7893 8010 8082 17503 96295 29477 13385 160658 4881 5546 5772 7623 11193 9169 9077 9787 10881 11887 10227 9264 8809 9127 9618 16199 97917 40492 21747
2 041010010 2 仙台市青葉区 青葉町 0 <NA> <NA> 649 16 15 17 23 53 62 49 40 40 40 45 33 28 38 38 48 413 143 67 307 8 10 7 10 26 30 27 21 18 14 25 16 11 15 16 25 198 60 29 342 8 5 10 13 27 32 22 19 22 26 20 17 17 23 22 23 215 83 38
3 041010020 2 仙台市青葉区 あけぼの町 0 <NA> <NA> 741 23 18 13 26 32 55 48 42 60 51 48 47 38 55 43 54 447 209 111 365 12 10 4 15 16 22 28 19 39 27 25 16 23 24 21 26 230 89 44 376 11 8 9 11 16 33 20 23 21 24 23 31 15 31 22 28 217 120 67
4 041010030 3 仙台市青葉区 旭ケ丘 0 <NA> <NA> 9160 279 289 272 315 766 880 771 643 633 713 561 493 436 363 358 840 6211 1671 950 4274 149 161 141 155 315 366 352 308 296 350 278 237 220 165 167 451 2877 697 365 4886 130 128 131 160 451 514 419 335 337 363 283 256 216 198 191 389 3334 974 585
--
5940 04606004015 4 南三陸町 歌津字石浜 0 <NA> <NA> 295 9 7 11 6 6 12 6 18 16 19 21 36 18 30 26 27 158 110 54 146 5 5 5 3 5 8 2 6 8 10 14 16 8 11 20 15 80 51 20 149 4 2 6 3 1 4 4 12 8 9 7 20 10 19 6 12 78 59 34
5941 04606004016 4 南三陸町 歌津字田の浦 0 <NA> <NA> 144 5 2 5 3 5 7 5 7 6 5 15 12 17 14 4 12 82 50 32 66 - 1 3 1 3 4 3 4 3 2 6 5 9 9 4 4 40 22 9 78 5 1 2 2 2 3 2 3 3 3 9 7 8 5 - 8 42 28 23
5942 04606004017 4 南三陸町 歌津字草木沢 0 <NA> <NA> 457 21 16 18 16 12 13 20 18 18 41 34 28 37 46 43 55 237 165 76 234 11 6 8 10 11 5 9 10 11 22 18 15 17 19 29 25 128 81 33 223 10 10 10 6 1 8 11 8 7 19 16 13 20 27 14 30 109 84 43
5943 04606004018 4 南三陸町 歌津字伊里前 0 <NA> <NA> - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
数値のはずの列に混じる - と X を NA 扱い
これでようやくスタート地点。。。
dfL = readr::read_csv(infile, locale = sjis, col_select = seq(1, 7)) |>
dplyr::slice(-1)
dfR = readr::read_csv(infile, locale = sjis, col_select = -seq(1, 7),
skip = 1L, na = c("-", "X"))
raw_miyagi = dplyr::bind_cols(dfL, dfR) |> print()
KEY_CODE HYOSYO CITYNAME NAME HTKSYORI HTKSAKI GASSAN 総数、年齢「不詳」含む 総数０〜４歳 総数５〜９歳 総数１０〜１４歳 総数１５〜１９歳 総数２０〜２４歳 総数２５〜２９歳 総数３０〜３４歳 総数３５〜３９歳 総数４０〜４４歳 総数４５〜４９歳 総数５０〜５４歳 総数５５〜５９歳 総数６０〜６４歳 総数６５〜６９歳 総数７０〜７４歳 総数１５歳未満 総数１５〜６４歳 総数６５歳以上 総数７５歳以上 男の総数、年齢「不詳」含む 男０〜４歳 男５〜９歳 男１０〜１４歳 男１５〜１９歳 男２０〜２４歳 男２５〜２９歳 男３０〜３４歳 男３５〜３９歳 男４０〜４４歳 男４５〜４９歳 男５０〜５４歳 男５５〜５９歳 男６０〜６４歳 男６５〜６９歳 男７０〜７４歳 男１５歳未満 男１５〜６４歳 男６５歳以上 男７５歳以上 女の総数、年齢「不詳」含む 女０〜４歳 女５〜９歳 女１０〜１４歳 女１５〜１９歳 女２０〜２４歳 女２５〜２９歳 女３０〜３４歳 女３５〜３９歳 女４０〜４４歳 女４５〜４９歳 女５０〜５４歳 女５５〜５９歳 女６０〜６４歳 女６５〜６９歳 女７０〜７４歳 女１５歳未満 女１５〜６４歳 女６５歳以上 女７５歳以上
1 04101 1 仙台市青葉区 <NA> 0 <NA> <NA> 311590 10231 11633 11838 15944 23772 17838 17677 19028 21113 23299 20544 18295 16702 17137 17700 33702 194212 69969 35132 150932 5350 6087 6066 8321 12579 8669 8600 9241 10232 11412 10317 9031 7893 8010 8082 17503 96295 29477 13385 160658 4881 5546 5772 7623 11193 9169 9077 9787 10881 11887 10227 9264 8809 9127 9618 16199 97917 40492 21747
2 041010010 2 仙台市青葉区 青葉町 0 <NA> <NA> 649 16 15 17 23 53 62 49 40 40 40 45 33 28 38 38 48 413 143 67 307 8 10 7 10 26 30 27 21 18 14 25 16 11 15 16 25 198 60 29 342 8 5 10 13 27 32 22 19 22 26 20 17 17 23 22 23 215 83 38
3 041010020 2 仙台市青葉区 あけぼの町 0 <NA> <NA> 741 23 18 13 26 32 55 48 42 60 51 48 47 38 55 43 54 447 209 111 365 12 10 4 15 16 22 28 19 39 27 25 16 23 24 21 26 230 89 44 376 11 8 9 11 16 33 20 23 21 24 23 31 15 31 22 28 217 120 67
4 041010030 3 仙台市青葉区 旭ケ丘 0 <NA> <NA> 9160 279 289 272 315 766 880 771 643 633 713 561 493 436 363 358 840 6211 1671 950 4274 149 161 141 155 315 366 352 308 296 350 278 237 220 165 167 451 2877 697 365 4886 130 128 131 160 451 514 419 335 337 363 283 256 216 198 191 389 3334 974 585
--
5940 04606004015 4 南三陸町 歌津字石浜 0 <NA> <NA> 295 9 7 11 6 6 12 6 18 16 19 21 36 18 30 26 27 158 110 54 146 5 5 5 3 5 8 2 6 8 10 14 16 8 11 20 15 80 51 20 149 4 2 6 3 1 4 4 12 8 9 7 20 10 19 6 12 78 59 34
5941 04606004016 4 南三陸町 歌津字田の浦 0 <NA> <NA> 144 5 2 5 3 5 7 5 7 6 5 15 12 17 14 4 12 82 50 32 66 NA 1 3 1 3 4 3 4 3 2 6 5 9 9 4 4 40 22 9 78 5 1 2 2 2 3 2 3 3 3 9 7 8 5 NA 8 42 28 23
5942 04606004017 4 南三陸町 歌津字草木沢 0 <NA> <NA> 457 21 16 18 16 12 13 20 18 18 41 34 28 37 46 43 55 237 165 76 234 11 6 8 10 11 5 9 10 11 22 18 15 17 19 29 25 128 81 33 223 10 10 10 6 1 8 11 8 7 19 16 13 20 27 14 30 109 84 43
5943 04606004018 4 南三陸町 歌津字伊里前 0 <NA> <NA> NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA
あとはdplyrとtidyrでお掃除
まだ罠はたくさん: 列名の頭にスペース。全角数字。変な区分。
tidy_miyagi = raw_miyagi |>
dplyr::filter(HYOSYO == 1) |>
dplyr::select(CITYNAME, matches("[男女].+[歳上]$")) |>
tidyr::pivot_longer(!CITYNAME, names_to = "category", values_to = "count") |>
dplyr::mutate(category = stringr::str_trim(category)) |>
tidyr::separate(category, c("sex", "age"), 1) |>
dplyr::mutate(age = stringi::stri_trans_nfkc(age)) |>
tidyr::separate(age, c("lower", "upper"), "〜", fill = "right") |>
dplyr::mutate(lower = readr::parse_number(lower),
upper = readr::parse_number(upper)) |>
dplyr::filter((upper - lower) < 5 | lower == 75) |>
dplyr::mutate(age = (lower + upper) / 2) |>
print()
CITYNAME sex lower upper count age
1 仙台市青葉区 男 0 4 5350 2
2 仙台市青葉区 男 5 9 6087 7
3 仙台市青葉区 男 10 14 6066 12
4 仙台市青葉区 男 15 19 8321 17
--
1245 南三陸町 女 60 64 507 62
1246 南三陸町 女 65 69 553 67
1247 南三陸町 女 70 74 450 72
1248 南三陸町 女 75 NA 1602 NA
ここまでくれば作図は簡単
tidy_miyagi |>
dplyr::mutate(count = ifelse(sex == "男", -1, 1) * count) |>
ggplot() +
geom_col(aes(age, count, fill = sex)) +
facet_wrap(vars(CITYNAME), nrow = 4L) +
coord_flip() + theme_minimal(base_size = 15)
前処理の前処理が必要になるようなデータは困る
Rを使ったデータ前処理を習得した✨
雑然データからでも下ごしらえして解析できるぞ！
とはいえ…
- まず、読み込めないと始まらない。
- dplyrやtidyrで手に負えないほど不規則では困る。
自分が一次データを作る立場なら、どう気をつけたらいいか？
総務省が2020年に発表して「全国民が見るべき」と話題になった
「統計表における機械判読可能なデータの表記方法の統一ルール」
から抜粋。
データ入力の心得 1. 1つのセルに1つの値
悪い形に遭遇したときに使える関数:
tidyr::separate(), tidyr::separate_rows(), stringr::str_split(), stringr::str_extract()
データ入力の心得 1. 1つのセルに1つの値
悪い形に遭遇したときに使える関数:
tidyr::separate(), tidyr::separate_rows(), stringr::str_split(), stringr::str_extract()
データ入力の心得 2. 余計な文字列が数値扱いを妨げる
悪い形に遭遇したときに使える関数:
readr::parse_number(), stringr::str_remove(), stringr::str_replace()
データ入力の心得 2. 余計な文字列が数値扱いを妨げる
悪い形に遭遇したときに使える関数:
readr::parse_number(), stringr::str_remove(), stringr::str_replace()
データ入力の心得 2. 余計な文字列が数値扱いを妨げる
いやいや、欄外に記載するのも不規則でダメですよ。別シートか、別カラムに。
悪い形に遭遇したときに使える関数:
readr::parse_number(), stringr::str_remove(), stringr::str_replace()
データ入力の心得 3. 空白や改行で見た目を整えない
機械にとって " A" と "A" は違う値。
悪い形に遭遇したときに使える関数:
stringr::str_remove(), stringr::str_replace()
データ入力の心得 3. 空白や改行で見た目を整えない
悪い形に遭遇したときに使える関数:
stringr::str_remove(), stringr::str_replace()
データ入力の心得 4. セルを結合しない、絶対に
悪い形に遭遇したときに使える関数:
tidyr::fill()
データ入力の心得 4. セルを結合しない、絶対に
悪い形に遭遇したときに使える関数:
tidyr::fill()
データ入力の心得 5. 同じものが続いても省略しない
人が見れば分かるけど、機械も同じように読んではくれない。
悪い形に遭遇したときに使える関数:
tidyr::fill(), tidyr::separate(), stringr::str_replace()
データ入力の心得 6. 機種依存文字を避ける
できれば全角・日本語も避けてすべて半角英数で。
悪い形に遭遇したときに使える関数:
stringi::stri_trans_nfkc()
データ入力の心得 7. データを分断しない
データ入力の心得 7. データを分断しない
データ入力の心得 8. 1つのシートには1つの表
もっと言うと、複数シートのファイルを作らず、それぞれ別ファイルに。
データ入力の心得 9. その他
- 1行目ではすべての列に名前を与える。
- なんとなく1行目は空行、1列目は名前無しとかやりがち。ダメ。
- 2行目以降に列名がハミ出してもダメ。
- プログラム内で参照される変数名として成立するように。
- アンダースコア_以外の記号・空白を入れないのが安全。
- 先頭は数字も避ける。
- 保存形式はプレーンテキスト推奨。Excel形式(.xlsx)は扱いにくい。
- タブ区切り(.tsv) or カンマ区切り(.csv)
- ファイル名、フォルダ名も機械処理するつもりで規則的に。
データ入力の心得おさらい
つまり、Rで前処理・可視化・解析するぞという気持ち:
- 1つのセルに1つの値
- 余計な文字列が数値扱いを妨げる
- 空白や改行で見た目を整えない
- セルを結合しない、絶対に
- 同じものが続いても省略しない
- 機種依存文字、というか半角英数以外を避ける
- データを分断しない
- 1つのシート(ファイル)には1つの表
🔰 総務省ルールの「修正後」にまだ残る改善点をなるべく多く列挙しよう。
Excelの親切あるいは裏切りに要注意
表形式のファイルを閲覧・作成するソフトとしてよく普及。
しかし親切設計がユーザーに牙を剥く！
- 文字列が勝手に日付になる:
22-4,4-14 - 遺伝子名が勝手に日付になる:
MARCH1,SEPT1
(2020年、これを避けるために改名された)
🔰 CSVファイルを読み書きして恐怖を体感しよう:
gene,label
MARCH1,22-4
SEPT1,4-14
- 上の内容のプレーンテキストを作って
excel.csvと名付ける - Excelでそれで開き、日付になることを確認。
- それを別名のCSVとして保存し、テキストの中身を確認してレポート。
今回のお品書き
✅ データ入力の心得
- 野良データを読み込む難しさ
- 自分が一次データを作るときに気をつけること
⬜ データ解釈の心得 (理論や手法の詳細には触れない)
- 誤差 (バイアス、ばらつき)
- 統計的仮説検定
- 因果関係と相関関係
Garbage in, garbage out
どんなに素晴らしい統計処理をしようとも、データがゴミなら結論もゴミ。
- 目的に応じてちゃんとデータを採ることが重要！
- 何がデータを悪くしてしまうのか、見ていこう
観測には必ず何らかの歪みがある
- 現象・母集団の性質を知りたいけど、そのものは測れない。
- 質的な問題 → 測れる側面だけ測るしかない
- 量的な問題 → 一部の標本を抽出して測るしかない
2種類の誤差
- 系統的誤差 systematic error / バイアス bias
- 一定の傾向をもって生じる誤差。原因を特定して対処。
- 偶然誤差 random error
- 偶然生じてしまうばらつき。計測を繰り返して対処。
e.g., 服を着たまま体重測定を10回繰り返す → どれも違う値
選択バイアス
データの採れかたが結果に傾向をもたらしてしまう。
- 生存(者)バイアス
- 生き残ったものだけがデータとして観測される。
- 帰還した多数の戦闘機✈️✈️の損傷の分布(右図)。
被弾 → 生還 → 観測される
被弾 → 墜落 → 観測されない
損傷が観察されなかった場所をむしろ補強すべき。 - サンプリングバイアス
- サンプル対象を決めた時点で偏っている。
- 📞電話でアンケート調査 → 電話機を持てる裕福な人、その時間に通話できる人、 知らない番号からでも出る人、調査している新聞社に共感している人。
🔰 生物学研究における選択バイアスにはどういうものがありそう？
測定基準に関するバイアスの例
- 基準が揃っていない
- 国民1人あたりの弁護士数を国際比較すると日本が少ない。
諸外国では税理士や司法書士を弁護士に含めて集計している。 - 基準が時間とともに変化
- 自閉症と診断された児童の数は見かけ上年々増えている。
診断基準の変更、病気の知名度向上が主な原因。
(親の高齢化の寄与を示唆する研究もあるらしい)
データの扱いに起因するバイアス
- チェリー・ピッキング 🍒
- 自分の仮説と矛盾する情報を無視し、都合のいい情報だけ集めること。
- 研究活動における3つの不正のうち改竄にあたる。
- 確証バイアス
- 無意識にチェリー・ピッキングしてしまう心理的傾向。
- 悪気なくやってしまわないよう、自覚を持つことが大事。
表に数字、裏にアルファベットが書いてあるカード。
「偶数の裏は必ずA、という法則があるか？」
これを確かめるには次の4枚のうちどれを裏返せばいいか？
たまたま生じたパターンにも意味を見出しがち
- 前後即因果の誤謬
- たまたま前後で生じた現象に因果関係があると思いがち
- ギャンブラーの誤謬
- N回も連続で裏が出たら、次は表が出やすいと思いがち
データ解釈に関わる認知バイアス
- 利用可能性ヒューリスティック
- 思い浮かびやすいものほど過大評価しがち
- e.g., 1文字目がkの単語、3文字目がkの単語、多いのはどっち？
- 代表性ヒューリスティック
- ステレオタイプを過大評価しがち
- e.g., リンダは31歳の独身女性。率直で聡明。大学では哲学専攻。
人種差別や社会正義などに深く関心を持ち、学生時代には反核デモにも参加。
現在のリンダは次のどちらである可能性が高いか？- リンダは銀行員として働いている。
- リンダは銀行員として働きながら、フェミニスト運動に参加している。
認知バイアス・誤謬はほかにもいろいろある。
データ解析とは関係なく、知っておくのは悪くない。
偶然誤差によるばらつき方 = 確率分布
次回 7. 統計モデリング1 でやります。
背後にあるメカニズムによっていろんな形になる。
分布を特徴づける代表値 central tendency
- 平均値 mean
- 和を観察数で割る
- 中央値 median
- 順に並べて真ん中
- 最頻値 mode
- 最も頻度が高い値
目的や状況に応じて使い分けよう。
- 外れ値に対する応答
- もし総資産額20兆円の大富豪が鳥取県に引っ越してきたら
→ 県民の平均資産は4000万円上昇。中央値・最頻値はほぼそのまま。
ばらつきを捉える記述統計量
- 分散 variance
- 平均値からの差の自乗の平均。 $\frac 1 n \sum _i ^n (X_i - \bar X)^2$
- これの平方根が標準偏差 (standard deviation)。
- Percentile, Quantile (四分位)
- 小さい順にならべて上位何%にあるか。
- 中央値 = 50th percentile = 第二四分位(Q2)
記述統計量に頼りすぎず分布を可視化する
同じデータでも見せ方で印象・情報量が変わる。
2つの量の関係性: 大小の比較
ばらつきの度合いも加味して判断する。
観測値1つずつ。
たまたまかも。
たまたまかも。
ばらつき大きい。
Bが高いのもたまたま?
Bが高いのもたまたま?
ばらつき小さい。
AとBには差がありそう。
AとBには差がありそう。
「こんなことがたまたま起こる確率はすごく低いです！」
をちゃんと示す手続きが統計的仮説検定→
目的やデータに応じて正しい検定を選ぶ
例えば「2つの量の大小を比較したい」だけでもいくつかある。
今回は選定には深入りせず、検定の考え方だけ触れておく。
統計的仮説検定の考え方
🎲サイコロを10回振ったら9回も6の目が出たぞ。イカサマじゃないか？
- 帰無仮説
- 6の目が出る確率 = 1/6。普通のサイコロです。
- 対立仮説
- 6の目が出る確率 ≠ 1/6。イカサマです。
\[\begin{split}
p = \binom {10} {9} \times {\frac 1 6} ^ {9} \times {\frac 5 6} ^ 1 + {\frac 1 6} ^ {10}
= 8.43 \times 10 ^ {-7}
\end{split}\]
- 帰無仮説のもとで今回のデータ以上に極端な値が得られる確率pを計算
- この p-value が有意水準(危険率)$\alpha$より低い場合、帰無仮説を棄却
(大概 $\alpha$ = 0.05 とか 0.01 とか。予め決めておく。) - 対立仮説を採択。6の目の出る確率は1/6と有意に異なる。
統計的仮説検定の考え方
🎲サイコロを12回振ったら4回も6の目が出たぞ。
これも1/6よりは高いけど、イカサマのダイスかな？
帰無仮説のもとで12回中4回以上6が出る確率: $p = 0.125 > \alpha$
帰無仮説を棄却できず。6の目の出る確率が1/6と有意に異なるとは言えない。
平たく言えば「普通のサイコロでもそんなに珍しい出目じゃない」
多重検定で偽陽性のリスクが上昇
危険率$\alpha=0.05$の検定は、最大5%の確率で誤って帰無仮説を棄却(偽陽性)。
同様の検定を10回やると、そのうち少なくとも1つで偽陽性になる確率は
$1 - (1 - 0.05)^{10} \approx 0.40$ まで上昇。
- 多重比較補正
- 検定を繰り返したぶんだけ基準を厳しくして偽陽性率を抑えること。
- Bonferroni法、Holm法、Benjamini and Hochberg法など。
広くて深い統計解析の世界
データの荒波を泳ぎきってもどこにも「究極の真実」などありはしないのだ。 統計学はそのときその場かぎりでの「最良の結論」を導くための方便に過ぎないのだ。 —— 三中信宏「統計思考の世界」
2つの量の関係性: 相関関係と因果関係
- 相関関係
- ある値が大きいほど、別の値も大きいor小さい。
- e.g., 数学の成績と物理の成績。
- 因果関係
- ある事象が別の事象に影響を与える。
- e.g., 1時間勉強するごとに成績が3点伸びる。
- 因果関係があれば相関関係として表れる
- 相関関係があるからといって因果関係もあるとは限らない→
因果・相関を見誤るパターン1: 交絡因子
誤: 🍦アイスの売り上げが増えるほど🍺ビールの売り上げが増える。
正: 🍦アイスも🍺ビールも気温が高いほど売れる。
因果・相関を見誤るパターン2: 逆の因果関係
誤: 警察官が多いほど犯罪が増える。
正: 犯罪が多いから警察官が多く配備される。
因果・相関を見誤るパターン3: 選択バイアス
(x + y) が一定の幅に収まるようなペアだけを集めてしまうとか
因果・相関を見誤るパターン4: 外れ値・グループ構造
少数の外れ値やグループ構造によって相関係数rが跳ね上がりうる。
因果・相関を見誤るパターン5: 偶然
ニコラス・ケイジの映画出演が増えるほど溺死する人が増える？
このウェブサイトにはジョークとしてこうした例が多数集められている↑
変数間の関係性まとめ
- 相関関係があるからといって因果関係もあるとは限らない。
- 逆の因果関係や見せかけの相関にも要注意。
🔰 これら4つの関係性に該当する事例をそれぞれ1つ以上探してみよう。
🔰 3日目の課題2
書籍、ニュースサイト、ブログ記事などを検索し、
良くないデータ解析の例を探して改善案を挙げよう。
- データの集め方(サンプルの選び方)が偏っていて良くない例
- データの解釈(相関や因果)を誤っている例
- 図の描き方で嘘をついている例 (誇張、印象操作)
- 生データを公開してくれるのはいいけど形式が汚くて使いにくい例
それぞれの例について、下記の項目を報告してください:
- 出典。URLや書籍名。
- データの概要、書き手の主張。
- ここが良くないので、こう改善すべき。
今日の残り時間
- 班やTAに相談し、消化しきれなかった部分をなるべく解消する。
- まずは個人で課題1と2に取り組み、ある程度できたら班内で見せ合う。
- 班員のよかった部分を自分のレポートに取り込む。
- 最もよかった解析をTAに報告して解散。
- 残ってほかの課題に取り組んでもOK。
参考文献
- 「Rを用いたデータ解析の基礎と応用」 石川由希 2022 名古屋大学
- データ解析のための統計モデリング入門 久保拓弥 2012
- 分析者のためのデータ解釈学入門 江崎貴裕 2020
- データ分析のための数理モデル入門 江崎貴裕 2020
- 統計学を哲学する 大塚淳 2020
- 統計思考の世界 三中信宏 2018
- R Markdown