LEC06: グラフの実践2（ggplot2）改訂版（11/8版）
Kenri Kodaka
2025-11-6
[準備1]
データフレームをインポート
◯◯◯ 使用するライブラリ
ggpplot2ライブラリ（グラフ描画パッケージ）およびdplyr（データ操作）を使用するために、tidyverseパッケージを、また、複数グラフ配置のためにpatchworkライブラリを
install.packages("tidyverse")
install.packages("patchwork")library(ggplot2)
library(dplyr)
library(patchwork)◯◯◯ データフレームをインポート
今回の演習では、集団フリーシャッター実験のデータを扱います。
実験の内容については以下の映像を参考にしてください。 - [YOUTUBE] 集団フリーシャッター課題
まずは以下のCSVを読み込みましょう。
# csvは改訂版
url = "https://lab.kenrikodaka.com/_download/csv/freeshutter_230601m.csv"
dat = read.csv(url)最初の6行を確認します。
head(dat)## seat a1 b1 b2 b3 seat_ix seat_iy ## 1 192 2.2 1.5 2.7 5.6 12 16 ## 2 168 0.8 4.4 10.0 14.7 12 14 ## 3 59 0.6 0.6 0.8 0.9 11 5 ## 4 60 2.9 3.4 5.3 9.9 12 5 ## 5 89 6.2 0.3 10.2 18.1 5 8 ## 6 180 0.1 0.1 0.5 0.8 12 15
個々の列に対応する変数の内容は以下の通りです。
- seat（座席番号）
- a1（実験Aのシャッター時間）
- b1（実験Bのシャッター時間：1回目）
- b2（実験Bのシャッター時間：2回目）
- b3（実験Bのシャッター時間：3回目）
- seat_ix（座席の左右位置：1-12）
- seat_iy（座席の奥行き位置：1-22）
なおこのデータは、以下の二種類の集団フリーシャッター実験を含んでいます。
- 実験A（20秒間の間に１度だけシャッターを押す）
- 実験B（20秒間の間に最大で３回シャッターを押す）
# 全部で119人分のデータがあることが確認できます。
dim(dat)## [1] 119 7
◯◯◯ pivot_longerによるデータの整形
現在のデータフレームの仕様は、一行に一人の実験参加者が対応し、それぞれの行が、4つのシャッター時間の計測値(a1,b1,b2,b3)を持ちます。
これを、
pivot_longerを用いて、一行に一つの計測時間を対応させるかたちに変形させます。この場合、列の変数（a1,b1,b2,b3)を、新たな変数「group」の値として、もとの計測値を「time」に対応させます。
library(tidyr)
dat.fs = dat
dat.fs = dat.fs %>%
pivot_longer(
cols = c(a1, b1, b2, b3), # 値変換する変数を指定
names_to = "group", # colsを扱う新たな変数
values_to = "time" # もともとの値を収納する変数
)
dat.fs## # A tibble: 476 × 5 ## seat seat_ix seat_iy group time ## <int> <int> <int> <chr> <dbl> ## 1 192 12 16 a1 2.2 ## 2 192 12 16 b1 1.5 ## 3 192 12 16 b2 2.7 ## 4 192 12 16 b3 5.6 ## 5 168 12 14 a1 0.8 ## 6 168 12 14 b1 4.4 ## 7 168 12 14 b2 10 ## 8 168 12 14 b3 14.7 ## 9 59 11 5 a1 0.6 ## 10 59 11 5 b1 0.6 ## # ℹ 466 more rows
最後にgroupを、実験の種類「exp」と撮影順序「order」に分解します。
dat.fs = dat.fs %>%
mutate(exp = if_else(group=="a1","Single","Triple")) %>%
mutate(order = case_when(group=="a1" ~ "1st",
group=="b1" ~ "1st",
group=="b2" ~ "2nd",
group=="b3" ~ "3rd")) %>%
select(-group) #groupを削除
dat.fs## # A tibble: 476 × 6 ## seat seat_ix seat_iy time exp order ## <int> <int> <int> <dbl> <chr> <chr> ## 1 192 12 16 2.2 Single 1st ## 2 192 12 16 1.5 Triple 1st ## 3 192 12 16 2.7 Triple 2nd ## 4 192 12 16 5.6 Triple 3rd ## 5 168 12 14 0.8 Single 1st ## 6 168 12 14 4.4 Triple 1st ## 7 168 12 14 10 Triple 2nd ## 8 168 12 14 14.7 Triple 3rd ## 9 59 11 5 0.6 Single 1st ## 10 59 11 5 0.6 Triple 1st ## # ℹ 466 more rows
[Graph1_ShutterHistogram]
シャッターのヒストグラム
◯◯◯ Graph1A（Single実験）
# Single実験のデータフレームを取得
dat_si = dat.fs %>%
filter(exp == "Single")
# [ggplot]
## x軸をtimeに設定
gp_a = ggplot(dat_si,aes(x=time))
# [geom_histogram] ヒストグラムの作成
## 一つのビンの幅を0.1sec、塗りつぶしを白、枠線を黒
gp_a = gp_a + geom_histogram(binwidth=0.1, fill="white",colour="black")
# Y軸（Count）の最大を10とする
gp_a = gp_a + scale_y_continuous(limits=c(0,10),
breaks=2*(0:5))
# X軸の範囲と目盛の位置を決める
gp_a = gp_a + scale_x_continuous(limits=c(0,20),
breaks=5*(0:4))
# [labs]
## グラフのタイトル、XY軸のタイトル
gp_a = gp_a + labs(x = "Time[s]", y = "Count")gp_a + theme_minimal() +
ggtitle("Graph1A_ShutterHistogram")
◯◯◯ Graph1B（Triple実験）
# Triple実験のデータフレームを取得
dat_tri = dat.fs %>%
filter(exp == "Triple")
# [ggplot]
## x軸をtimeに設定
gp_b = ggplot(dat_tri,aes(x=time))
gp_b = gp_b + geom_histogram(binwidth=0.1, fill="white",colour="black")
gp_b = gp_b + scale_y_continuous(limits=c(0,10),breaks=2*(0:5))
gp_b = gp_b + scale_x_continuous(limits=c(0,20),breaks=5*(0:4))
gp_b = gp_b + labs(x = "Time[s]", y = "Count")gp_b + theme_minimal() +
ggtitle("Graph1B_ShutterHistogram")
◯◯◯ Graph1C（Triple実験、順序別、積み上げ方式）
# Triple実験のデータフレームを取得
dat_tri = dat.fs %>%
filter(exp == "Triple")
# [ggplot]
## x軸をtimeに設定, 塗りつぶしをorderに設定
gp_c = ggplot(dat_tri,aes(x=time,fill=order))
# 積み上げ順序の反転（factor順序：下から上へ）
## 1->2->3の順に積み上げるようにする
gp_c = gp_c + geom_histogram(binwidth=0.1,colour="black",
position = position_stack(reverse = TRUE))
gp_c = gp_c + scale_y_continuous(limits=c(0,10),breaks=2*(0:5))
gp_c = gp_c + scale_x_continuous(limits=c(0,20),breaks=5*(0:4))
# labs()で各種のタイトルをまとめて変更（fillは凡例）
gp_c = gp_c + labs(x = "Time[s]", y = "Count", fill="Order")
# 凡例の文字列（変更なし：実行しなくても良い）
gp_c = gp_c + scale_fill_discrete(
limits=c("1st","2nd","3rd"),labels=c("1st","2nd","3rd"))gp_c + theme_minimal() +
ggtitle("Graph1C_ShutterHistogram")
◯◯◯ Graph1D（Triple実験、順序別）
# Triple実験のデータフレームを取得
dat_tri = dat.fs %>%
filter(exp == "Triple")
gp_d = ggplot(dat_tri,aes(x=time,fill=order))
gp_d = gp_d + geom_histogram(binwidth=0.1,colour="black")
gp_d = gp_d + scale_x_continuous(limits=c(0,20),breaks=2*(0:10))
gp_d = gp_d + scale_y_continuous(limits=c(0,10),breaks=5*(0:4))
# 凡例の文字列（変更なし：実行しなくても良い）
gp_d = gp_d + scale_fill_discrete(
limits=c("1st","2nd","3rd"),labels=c("1st","2nd","3rd"))
# labs()で各種のタイトルをまとめて変更（fillは凡例）
gp_d = gp_d + labs(x="Time[s]", y="Count", fill="Order")
# facetをgroupで分割、factの背景なども整える
gp_d = gp_d + facet_grid(order ~ .)gp_d + theme_minimal() +
ggtitle("Graph1D_ShutterHistogram")
[Graph2_ShutterBox]
シャッター時間のボックスプロット
シャッター時間をボックスプロットで描画します。
まずは、一つのグラフで出力する方法から。
# [ggplot]
## x軸をexp、y軸にtime, 塗りつぶしをtroupに設定
gp_a = ggplot(dat.fs,aes(x=exp,y=time,fill=order))
# [geom_boxplot]ボックスプロットによる描画
gp_a = gp_a + geom_boxplot()
gp_a = gp_a + scale_x_discrete(limits=c("Single","Triple"),
labels=c("Type Single", "Type Triple"))
gp_a = gp_a + scale_y_continuous(limits=c(0,20),breaks=5*(0:4))
# Bluesパレットを使って、fillを描画。
gp_a = gp_a + scale_fill_brewer(palette = "Blues")
# RColorBrewer::display.brewer.all()
gp_a = gp_a + labs(x = "Experiment", y = "Time[s]",fill="Order")gp_a + theme_minimal() +
ggtitle("Graph2_ShutterBox")
別の方法として、facetを使ってグラフをExpのタイプごとに切り分けます。
# [ggplot]
## x軸をorder、y軸にtime, 塗りつぶしをorderに設定
gp_b = ggplot(dat.fs,aes(x=order,y=time,fill=order))
# [geom_boxplot]ボックスプロットによる描画
gp_b = gp_b + geom_boxplot()
gp_b = gp_b + scale_x_discrete()
gp_b = gp_b + scale_y_continuous(limits=c(0,20),breaks=5*(0:4))
# Bluesパレットを使って、fillを描画。
gp_b = gp_b + scale_fill_brewer(palette = "Blues")
gp_b = gp_b + labs(x = "Order", y = "Time[s]",fill="Order")
# facetをgroupで分割、factの背景なども整える
gp_b = gp_b + facet_grid(. ~ exp)gp_b + theme_minimal() +
ggtitle("Graph2B_ShutterBox")
[Graph3_ShutterComma]
コンマ何秒でシャッターを押していたか
あらためて、dat.fsの構造は以下の通り。
head(dat.fs); dim(dat.fs)## # A tibble: 6 × 6 ## seat seat_ix seat_iy time exp order ## <int> <int> <int> <dbl> <chr> <chr> ## 1 192 12 16 2.2 Single 1st ## 2 192 12 16 1.5 Triple 1st ## 3 192 12 16 2.7 Triple 2nd ## 4 192 12 16 5.6 Triple 3rd ## 5 168 12 14 0.8 Single 1st ## 6 168 12 14 4.4 Triple 1st
## [1] 476 6
◯◯◯ 欠損値の削除
is.na()を使って、計測値が欠損している行を取り出してみます。 その結果、476行の中には42の欠損データがあり、 その多くが、triple実験の2~3回目であることがわかります。
dat.fs %>% filter(is.na(time))## # A tibble: 42 × 6 ## seat seat_ix seat_iy time exp order ## <int> <int> <int> <dbl> <chr> <chr> ## 1 156 12 13 NA Triple 3rd ## 2 204 12 17 NA Triple 2nd ## 3 204 12 17 NA Triple 3rd ## 4 227 11 19 NA Triple 2nd ## 5 227 11 19 NA Triple 3rd ## 6 37 1 4 NA Triple 3rd ## 7 38 2 4 NA Triple 2nd ## 8 38 2 4 NA Triple 3rd ## 9 94 10 8 NA Triple 3rd ## 10 236 8 20 NA Triple 2nd ## # ℹ 32 more rows
sum(is.na(dat.fs$time))## [1] 42
今後のために、dat.fsから欠損値を除いておきます（434行）。
dat.fs = dat.fs %>%
filter(!is.na(time))
dim(dat.fs)## [1] 434 6
◯◯◯ コンマN秒ごとに集計
dat.fs2 = dat.fs %>%
mutate(time2 = as.integer(floor(time*10)%%10)) %>%
select(time,time2,exp,order)
dat.fs2## # A tibble: 434 × 4 ## time time2 exp order ## <dbl> <int> <chr> <chr> ## 1 2.2 2 Single 1st ## 2 1.5 5 Triple 1st ## 3 2.7 7 Triple 2nd ## 4 5.6 6 Triple 3rd ## 5 0.8 8 Single 1st ## 6 4.4 4 Triple 1st ## 7 10 0 Triple 2nd ## 8 14.7 7 Triple 3rd ## 9 0.6 6 Single 1st ## 10 0.6 6 Triple 1st ## # ℹ 424 more rows
dat_comma = dat.fs2 %>%
group_by(exp,time2) %>% #実験の種類、コンマn秒でグループ分け
summarize(count=n()) #グループに基づき総数を数え上げる
dat_comma## # A tibble: 20 × 3 ## # Groups: exp [2] ## exp time2 count ## <chr> <int> <int> ## 1 Single 0 17 ## 2 Single 1 5 ## 3 Single 2 16 ## 4 Single 3 13 ## 5 Single 4 10 ## 6 Single 5 11 ## 7 Single 6 11 ## 8 Single 7 12 ## 9 Single 8 14 ## 10 Single 9 10 ## 11 Triple 0 42 ## 12 Triple 1 25 ## 13 Triple 2 32 ## 14 Triple 3 32 ## 15 Triple 4 35 ## 16 Triple 5 28 ## 17 Triple 6 20 ## 18 Triple 7 31 ## 19 Triple 8 33 ## 20 Triple 9 37
◯◯◯ Graph3A（X軸を0:9の順で）
gp_a = ggplot(dat_comma,aes(x=time2,y=count,colour=exp))
gp_a = gp_a + geom_line(size=0.5) + geom_point(size=3)
gp_a = gp_a + scale_colour_brewer(palette ="Set2")
gp_a = gp_a + scale_x_continuous(limits=c(-0.5,9.5),breaks=0:9)
gp_a = gp_a + scale_y_continuous(limits=c(0,50),breaks=seq(0,50,by=10))
gp_a = gp_a + labs(x="Comma.N",y="Count",colour="Exp Type")gp_a + theme_minimal() +
ggtitle("Graph3A_ShutterComma")
◯◯◯ Graph3B（X軸を6,7,8,9,0,1,2,3,4,5の順で）
# time2を指定した順序でfactor化
order_n = c(6,7,8,9,0,1,2,3,4,5)
dat_comma = dat_comma %>%
mutate(time2 = factor(time2, levels = order_n))# groupで線を引くグループを指定します。
gp_b = ggplot(dat_comma,
aes(x=time2,y=count,group=exp,colour=exp))
## これは、後続するgeom_line()が
##「同じグループ内の点を順に線でつなぐ」構造になっているためです
# x=5（time2=0に相当）で縦の点線を引きます。
gp_b = gp_b + geom_vline(xintercept = 5, linetype = "dotted",
size=0.5,colour="black")
gp_b = gp_b + scale_colour_brewer(palette ="Set2")
gp_b = gp_b + geom_line(size=0.5) + geom_point(size=3)
# xはfactorのため、離散値（discrete）として扱われることに注意。
## paste0は、文字列の接合をします。
gp_b = gp_b + scale_x_discrete(labels=paste0(".",c(6:9,0:5)))
gp_b = gp_b + scale_y_continuous(limits=c(0,50),breaks=seq(0,50,by=10))
gp_b = gp_b + labs(x="Comma.N",y="Count",colour="Exp Type");gp_b + theme_minimal() +
ggtitle("Graph3B_ShutterComma")
[Graph4_Correlation]
シャッターとシャッターの時間差の相関（Triple実験）
まずは、
pivot_longerをする前のデータを確認します。
head(dat,10)## seat a1 b1 b2 b3 seat_ix seat_iy ## 1 192 2.2 1.5 2.7 5.6 12 16 ## 2 168 0.8 4.4 10.0 14.7 12 14 ## 3 59 0.6 0.6 0.8 0.9 11 5 ## 4 60 2.9 3.4 5.3 9.9 12 5 ## 5 89 6.2 0.3 10.2 18.1 5 8 ## 6 180 0.1 0.1 0.5 0.8 12 15 ## 7 223 14.8 0.8 9.8 19.5 7 19 ## 8 213 0.7 1.3 2.4 3.9 9 18 ## 9 214 1.9 1.2 2.2 3.1 10 18 ## 10 191 0.3 18.2 18.5 18.9 11 16
オリジナルの構造を保ったまま、 Triple実験の時間差のみを抽出するデータフレーム（
dat_dt）を作成します。
# datをdat_dtにコピーします。
dat_dt = dat
dat_dt = dat_dt %>%
mutate(time1 = b1, #変数名をtime1にコピー
time2 = b2, #同様
time3 = b3, #同様
time21 = time2-time1, #シャッター1とシャッター2の時間差
time32 = time3-time2 #シャッター2とシャッター3の時間差
) %>%
select(time1,time21,time32) #最初のシャッター時間、および時間差のみに注目
head(dat_dt)## time1 time21 time32 ## 1 1.5 1.2 2.9 ## 2 4.4 5.6 4.7 ## 3 0.6 0.2 0.1 ## 4 3.4 1.9 4.6 ## 5 0.3 9.9 7.9 ## 6 0.1 0.4 0.3
◯◯◯ Graph4A（time1 と time21の相関図）
まずは、time1とtime21の相関をグラフ化してみます。
gp_a = ggplot(dat_dt, aes(x = time1, y = time21))
# [geom_point] 散布図（同一行のxとyの関係）
gp_a = gp_a + geom_point(size=2,shape=21)
# [geom_smooth] 線形回帰直線
gp_a = gp_a + geom_smooth(method="lm",fomula='y~x')
gp_a = gp_a + scale_x_continuous(limits=c(0,10),breaks=c(0,5,10))
gp_a = gp_a + scale_y_continuous(limits=c(0,10),breaks=c(0,5,10))
# [geom_segment] 任意の直線
# (x,y) to (xend,yend)の座標間に線を引く
gp_a = gp_a + geom_segment(x=0,y=0,xend=10,yend=10,
size=0.5,linetype="dotted")
gp_a = gp_a + labs(x = "Time1", y = "Time1 to Time2")最初にシャッターを押すまでの待機時間は、次にシャッターを押すまでの時間と まるで相関していないことがわかります。
gp_a = gp_a + theme_minimal() +
ggtitle("Graph4A_Correlation1")
gp_a
◯◯◯ Graph4B（time21 と time32の相関図）
time21とtime32の相関も確認してみましょう。
gp_b = ggplot(dat_dt, aes(x = time21, y = time32))
gp_b = gp_b + geom_point(size=2,shape=21)
gp_b = gp_b + geom_smooth(method="lm",fomula='y~x')
gp_b = gp_b + scale_x_continuous(limits=c(0,10),breaks=c(0,5,10))
gp_b = gp_b + scale_y_continuous(limits=c(0,10),breaks=c(0,5,10))
gp_b = gp_b + geom_segment(x=0,y=0,xend=10,yend=10,
size=0.5,linetype="dotted")
gp_b = gp_b + labs(x = "Time1 to Time2", y = "Time2 to Time3")シャッター間の時間差にはかなり強い相関が見られることがわかります。
gp_b = gp_b + theme_minimal() +
ggtitle("Graph4B_Correlation2")
gp_a | gp_b
（発展として）、順位相関係数を以下で求めることができます。
cor.test(dat_dt$time1, dat_dt$time21, method = "spearman")## ## Spearman's rank correlation rho ## ## data: dat_dt$time1 and dat_dt$time21 ## S = 208358, p-value = 0.6126 ## alternative hypothesis: true rho is not equal to 0 ## sample estimates: ## rho ## -0.04974421
cor.test(dat_dt$time21, dat_dt$time32, method = "spearman")## ## Spearman's rank correlation rho ## ## data: dat_dt$time21 and dat_dt$time32 ## S = 71330, p-value = 1.905e-05 ## alternative hypothesis: true rho is not equal to 0 ## sample estimates: ## rho ## 0.4319928
[Graph5_DirectionEffect]
隣接空間（方向）の影響
前の人に影響されるか？
後ろの人に影響されるか？
今回は、single実験のみを扱います。
まずは解析対象となる、座席（seatx,seaty）と計測値（time）のみのデータフレームを作ります。元となる
dat.fsの形式をまずは確認してください。
head(dat.fs)## # A tibble: 6 × 6 ## seat seat_ix seat_iy time exp order ## <int> <int> <int> <dbl> <chr> <chr> ## 1 192 12 16 2.2 Single 1st ## 2 192 12 16 1.5 Triple 1st ## 3 192 12 16 2.7 Triple 2nd ## 4 192 12 16 5.6 Triple 3rd ## 5 168 12 14 0.8 Single 1st ## 6 168 12 14 4.4 Triple 1st
filter→mutate→selectで、single実験の座席と時間のみを取り出し、これを新たに
dat_seatとします。
dat_seat = dat.fs %>%
filter(exp=="Single") %>%
filter(!is.na(seat)) %>% #seatが不明なものを取り除く。
mutate(seatx = seat_ix, #変数名を変えます。
seaty = seat_iy) %>%
select(seatx,seaty,time)
dat_seat## # A tibble: 117 × 3 ## seatx seaty time ## <int> <int> <dbl> ## 1 12 16 2.2 ## 2 12 14 0.8 ## 3 11 5 0.6 ## 4 12 5 2.9 ## 5 5 8 6.2 ## 6 12 15 0.1 ## 7 7 19 14.8 ## 8 9 18 0.7 ## 9 10 18 1.9 ## 10 11 16 0.3 ## # ℹ 107 more rows
まずは隣接関係にある者たちのシャッター時間を算出する関数を作成します。 引数は、自分の座席のxy座標です。
- timeForward：前方（最大）3人の平均シャッター時間の算出
- timeNeighbor：両隣（最大）2人の平均シャッター時間の算出
- timeBackward：後方（最大）3人の平均シャッター時間の算出
timeForward = function(sx,sy){
dat_seat %>%
filter(seatx>=sx-1 & seatx<=sx+1 & seaty==sy-1) %>%
summarize(avg = mean(time)) %>%
pull(avg)
}
timeNeighbor = function(sx,sy){
dat_seat %>%
filter((seatx==sx-1 | seatx==sx+1) & seaty==sy) %>%
summarize(avg = mean(time)) %>%
pull(avg)
}
timeBackward = function(sx,sy){
dat_seat %>%
filter(seatx>=sx-1 & seatx<=sx+1 & seaty==sy+1) %>%
summarize(avg = mean(time)) %>%
pull(avg)
}これらの関数を使って、各参加者に隣接する3領域のシャッター時間を取得します。この場合、以下のように計算できたらいいのですが、実際にはうまくいきません。それは、上で定義した関数の引数が、長さ１の値しか想定できていないからです。
dat_seat %>%
mutate(
time.f = timeForward(seatx,seaty),
time.n = timeNeighbor(seatx,seaty),
time.b = timeBackward(seatx,seaty)
)## # A tibble: 117 × 6 ## seatx seaty time time.f time.n time.b ## <int> <int> <dbl> <dbl> <dbl> <dbl> ## 1 12 16 2.2 NaN NaN NaN ## 2 12 14 0.8 NaN NaN NaN ## 3 11 5 0.6 NaN NaN NaN ## 4 12 5 2.9 NaN NaN NaN ## 5 5 8 6.2 NaN NaN NaN ## 6 12 15 0.1 NaN NaN NaN ## 7 7 19 14.8 NaN NaN NaN ## 8 9 18 0.7 NaN NaN NaN ## 9 10 18 1.9 NaN NaN NaN ## 10 11 16 0.3 NaN NaN NaN ## # ℹ 107 more rows
結論から言うと、この問題に対処するためには、
mapply関数を使う必要があります。以下が、その解決方法です。
dat_seat %>%
mutate(
time.f = mapply(timeForward,seatx,seaty),
time.n = mapply(timeNeighbor,seatx,seaty),
time.b = mapply(timeBackward,seatx,seaty)
)## # A tibble: 117 × 6 ## seatx seaty time time.f time.n time.b ## <int> <int> <dbl> <dbl> <dbl> <dbl> ## 1 12 16 2.2 0.1 0.3 2.2 ## 2 12 14 0.8 19.1 20 0.1 ## 3 11 5 0.6 2 2.9 NaN ## 4 12 5 2.9 2 0.6 NaN ## 5 5 8 6.2 11.6 10.3 10.7 ## 6 12 15 0.1 10.4 NaN 1.25 ## 7 7 19 14.8 5.4 16 0.4 ## 8 9 18 0.7 7.4 1.9 NaN ## 9 10 18 1.9 7.4 0.7 11 ## 10 11 16 0.3 5.2 2.2 2.2 ## # ℹ 107 more rows
ここで、一旦、
mapply関数について解説します。
mapply関数
mapply() は、複数のベクトルを並行して関数に適用し、それぞれの要素ペアに対して処理を繰り返すベースRのベクトル化関数です。ここでは典型的な使用例を確認します。
ifを含むようなある種の関数は、引数に（長さ１より大きな）ベクトルを持てません。 そのため、以下の計算はエラーとなります。
plus = function(x, y){
if(x+y<25) "A" else "B"
}
plus(1:10,11:20)## Error in if (x + y < 25) "A" else "B": the condition has length > 1mapply関数を使うと、長さ１の引数を前提とした関数に対しても、 繰り返し同じ関数を適用し、計算結果をベクトルとしてまとめられます。
# 以下の計算結果は
## c(plus(1,11),plus(2,12),plus(3,13),...)と同じ。
mapply(plus, 1:10, 11:20)## [1] "A" "A" "A" "A" "A" "A" "A" "B" "B" "B"
以上の結果は、
plus(x,y)をデータフレームの各行に対して適用した場合も同様です。 ここでは、seatxとseatyを加算したseatxyを新たに作成したいとします。mutateの中で、そのままplus関数を適用するとエラーになりますが、mapplyを使うことで、ベクトル化することができます。
dat_seat %>%
mutate(seatxy = plus(seatx,seaty))## Error in `mutate()`:
## ℹ In argument: `seatxy = plus(seatx, seaty)`.
## Caused by error in `if (x + y < 25) ...`:
## ! the condition has length > 1dat_seat %>%
mutate(seatxy = mapply(plus,seatx,seaty))## # A tibble: 117 × 4 ## seatx seaty time seatxy ## <int> <int> <dbl> <chr> ## 1 12 16 2.2 B ## 2 12 14 0.8 B ## 3 11 5 0.6 A ## 4 12 5 2.9 A ## 5 5 8 6.2 A ## 6 12 15 0.1 B ## 7 7 19 14.8 B ## 8 9 18 0.7 B ## 9 10 18 1.9 B ## 10 11 16 0.3 B ## # ℹ 107 more rows
mapply関数と同じことは、 rowwise関数で行単位でグループ化しても実現できます。
dat_seat %>%
rowwise() %>%
mutate(seatxy = plus(seatx,seaty)) ## # A tibble: 117 × 4 ## # Rowwise: ## seatx seaty time seatxy ## <int> <int> <dbl> <chr> ## 1 12 16 2.2 B ## 2 12 14 0.8 B ## 3 11 5 0.6 A ## 4 12 5 2.9 A ## 5 5 8 6.2 A ## 6 12 15 0.1 B ## 7 7 19 14.8 B ## 8 9 18 0.7 B ## 9 10 18 1.9 B ## 10 11 16 0.3 B ## # ℹ 107 more rows
このmapply関数を使うことで、簡単に既存のデータフレームを拡張して、 各席に対する、前方・近傍・後方の平均値を算出することができます。
dat_seat = dat_seat %>%
mutate(
time.f = mapply(timeForward,seatx,seaty),
time.n = mapply(timeNeighbor,seatx,seaty),
time.b = mapply(timeBackward,seatx,seaty)
)
dat_seat## # A tibble: 117 × 6 ## seatx seaty time time.f time.n time.b ## <int> <int> <dbl> <dbl> <dbl> <dbl> ## 1 12 16 2.2 0.1 0.3 2.2 ## 2 12 14 0.8 19.1 20 0.1 ## 3 11 5 0.6 2 2.9 NaN ## 4 12 5 2.9 2 0.6 NaN ## 5 5 8 6.2 11.6 10.3 10.7 ## 6 12 15 0.1 10.4 NaN 1.25 ## 7 7 19 14.8 5.4 16 0.4 ## 8 9 18 0.7 7.4 1.9 NaN ## 9 10 18 1.9 7.4 0.7 11 ## 10 11 16 0.3 5.2 2.2 2.2 ## # ℹ 107 more rows
◯◯◯ Graph5A（棒グラフ：LATE vs EARLY）
最後に3区間（前方・隣接・後方）の撮影時間と比べて 早くシャッターを押したか（EARLY）か遅く押したか（LATE）について、 mutate関数で、新たに属性を追加し並べていきます。
dat_seat = dat_seat %>%
mutate(order.f = if_else(time - time.f>=0,"LATE","EARLY"),
order.n = if_else(time - time.n>=0,"LATE","EARLY"),
order.b = if_else(time - time.b>=0,"LATE","EARLY"))
dat_seat[1:20,]## # A tibble: 20 × 9 ## seatx seaty time time.f time.n time.b order.f order.n order.b ## <int> <int> <dbl> <dbl> <dbl> <dbl> <chr> <chr> <chr> ## 1 12 16 2.2 0.1 0.3 2.2 LATE LATE LATE ## 2 12 14 0.8 19.1 20 0.1 EARLY EARLY LATE ## 3 11 5 0.6 2 2.9 NaN EARLY EARLY <NA> ## 4 12 5 2.9 2 0.6 NaN LATE LATE <NA> ## 5 5 8 6.2 11.6 10.3 10.7 EARLY EARLY EARLY ## 6 12 15 0.1 10.4 NaN 1.25 EARLY <NA> EARLY ## 7 7 19 14.8 5.4 16 0.4 LATE EARLY LATE ## 8 9 18 0.7 7.4 1.9 NaN EARLY EARLY <NA> ## 9 10 18 1.9 7.4 0.7 11 EARLY LATE EARLY ## 10 11 16 0.3 5.2 2.2 2.2 EARLY EARLY EARLY ## 11 9 17 7.4 5.2 NaN 1.3 LATE <NA> LATE ## 12 9 16 9.8 5.75 0.6 7.4 LATE LATE LATE ## 13 7 6 2.5 10 6.1 0.7 EARLY EARLY LATE ## 14 8 6 6.1 10 2.5 2.2 EARLY LATE LATE ## 15 10 7 3 NaN 5.1 0.8 <NA> EARLY LATE ## 16 9 7 3.7 6.1 1.85 1.73 EARLY LATE LATE ## 17 12 13 19.1 1.9 NaN 10.4 LATE <NA> LATE ## 18 12 17 2.2 1.25 NaN NaN LATE <NA> <NA> ## 19 9 11 9.7 3.2 16.8 9.3 LATE EARLY LATE ## 20 10 11 16.8 5.8 5.6 5.6 LATE LATE LATE
以下のように解析することができます。
- 席(5,8)は前方・両隣・後方よりも早くシャッターを押している。
- 席(7,10)は両隣よりも早くシャッターを押し、後方よりも遅くシャッターを押している。
# 前方との比較：ほぼ変わらない
table(dat_seat$order.f)## ## EARLY LATE ## 49 46
# 両隣との比較：ほぼ変わらない
table(dat_seat$order.n)## ## EARLY LATE ## 41 42
# 後方との比較：後方よりも早く押す傾向
table(dat_seat$order.b)## ## EARLY LATE ## 46 51
dat_dir = dat_seat %>%
pivot_longer(
cols = c(order.f,order.n,order.b),
names_to = "direction",
values_to = "dorder"
) %>%
filter(!is.na(dorder)) %>%
mutate(
direction = case_when(
direction=="order.f" ~ "FORWARD",
direction=="order.n" ~ "NEIGHBOR",
direction=="order.b" ~ "BACKWARD"
)
)
dat_dir## # A tibble: 275 × 8 ## seatx seaty time time.f time.n time.b direction dorder ## <int> <int> <dbl> <dbl> <dbl> <dbl> <chr> <chr> ## 1 12 16 2.2 0.1 0.3 2.2 FORWARD LATE ## 2 12 16 2.2 0.1 0.3 2.2 NEIGHBOR LATE ## 3 12 16 2.2 0.1 0.3 2.2 BACKWARD LATE ## 4 12 14 0.8 19.1 20 0.1 FORWARD EARLY ## 5 12 14 0.8 19.1 20 0.1 NEIGHBOR EARLY ## 6 12 14 0.8 19.1 20 0.1 BACKWARD LATE ## 7 11 5 0.6 2 2.9 NaN FORWARD EARLY ## 8 11 5 0.6 2 2.9 NaN NEIGHBOR EARLY ## 9 12 5 2.9 2 0.6 NaN FORWARD LATE ## 10 12 5 2.9 2 0.6 NaN NEIGHBOR LATE ## # ℹ 265 more rows
## x軸をdirectionに、塗りつぶしをdorderに設定
gp_a = ggplot(dat_dir,aes(x=direction, fill=dorder))
gp_a = gp_a + geom_bar(position=position_dodge(width=0.9), colour="black")
gp_a = gp_a + scale_x_discrete(limits=c("FORWARD","NEIGHBOR","BACKWARD"))
gp_a = gp_a + scale_y_continuous(limits=c(0,60))
gp_a = gp_a + scale_fill_brewer(palette="Set2", limits=c("EARLY","LATE"))
gp_a = gp_a + labs(x = "DIRECTION", y = "COUNT",fill="ORDER")gp_a + theme_minimal() +
ggtitle("Graph5A_DirectionEffect")
カイ二乗分布
統計で頻度の偏りの検定を行うのには、カイ二乗分布を使うことができます。 関数
chisq.testを使うと、以下のように一気に計算できます。
chisq.test(dat_dir$direction, dat_dir$dorder, correct=F) ## ## Pearson's Chi-squared test ## ## data: dat_dir$direction and dat_dir$dorder ## X-squared = 0.33183, df = 2, p-value = 0.8471
結果、
p>0.05で、方向による影響は統計的には認められないという結論となります。
◯◯◯ Graph5B（ボックスプロット）
新たに、隣接領域のシャッター時間と自分自身のシャッター時間の差分を
dif変数として列に追加します。difが正のとき「自分より遅いシャッター」、difが負のとき「自分より早くシャッター」を押していることになります。
dat_dir = dat_dir %>%
mutate(dif = case_when(
direction=="FORWARD" ~ time.f - time,
direction=="NEIGHBOR" ~ time.n - time,
direction=="BACKWARD" ~ time.b - time)
)
# グラフには使いませんが、平均値と中央値を比較します。
dat_dir %>%
group_by(direction) %>%
summarize(
avg = mean(dif), #平均値
med = median(dif) #中央値
)## # A tibble: 3 × 3 ## direction avg med ## <chr> <dbl> <dbl> ## 1 BACKWARD -0.124 0 ## 2 FORWARD 0.210 0.2 ## 3 NEIGHBOR 0.0711 -0.4
# [ggplot]
## x軸をdirection、y軸にdif, 塗りつぶしをdirectionに設定
gp_b = ggplot(dat_dir,aes(x=direction,y=dif,fill=direction))
# [geom_boxplot]ボックスプロットによる描画
gp_b = gp_b + geom_boxplot()
gp_b = gp_b + scale_x_discrete(limits=c("FORWARD","NEIGHBOR","BACKWARD"))
gp_b = gp_b + scale_y_continuous(limits=c(-6,6),breaks=2*(-3:3))
# Bluesパレットを使って、fillを描画。
gp_b = gp_b + scale_fill_brewer(palette = "Blues",
limits=c("FORWARD","NEIGHBOR","BACKWARD"))
gp_b = gp_b + labs(x = "DIRECTION", y = "Time[s]",fill="DIRECTION")gp_b + theme_minimal() +
ggtitle("Graph5B_DirectionEffect")
◯◯◯ Graph5C（ヒストグラム）
# 最初に、directionの並べ順を決めます（facet_gridの並べ順をFNBの順にします）
dat_dir$direction =
factor(dat_dir$direction,
levels = c("FORWARD","NEIGHBOR","BACKWARD"))
# [ggplot]
## x軸をdif、 塗りつぶしをdirectionに設定
gp_c = ggplot(dat_dir,aes(x=dif,fill=direction))
# x=0で縦の点線を引きます。
gp_c = gp_c + geom_vline(xintercept = 0, linetype = "solid",
size=0.4,colour="black")
# [geom_histogram]ヒストグラムによる描画
gp_c = gp_c + geom_histogram(binwidth=0.3,colour="black")
gp_c = gp_c + scale_x_continuous()
gp_c = gp_c + scale_y_continuous()
# Bluesパレットを使って、fillを描画。
gp_c = gp_c + scale_fill_brewer(palette = "Blues")
gp_c = gp_c + labs(x = "TIME DIFFERENCE", y = "COUNT",fill="DIRECTION")
gp_c = gp_c + facet_grid(direction ~ .)gp_c + theme_minimal() +
ggtitle("Graph5C_DirectionEffect")
[Graph6_NeighborEffect]
近傍性が同期率に与える影響
近くにいることで似たタイミングで押す傾向は強くなるか？
Graph6では、参加者間の近さ（近接性）が、シャッターを押すタイミングの近さ（同期性）に与える影響を検証します。素朴に考えると、参加者間の距離が近いほど、シャッタータイミングが近いと考えられます（以下ではこれを「近傍同期効果」と呼びます）。今回は、その仮説を検証してみましょう。
まずは、single実験の
seatxとseaty、timeを取り出します（seatが不明なものは除外）
dat_seat = dat.fs %>%
filter(exp=="Single") %>%
filter(!is.na(seat)) %>% #seatが不明なものを取り除く。
mutate(seatx = seat_ix, #変数名を変えます。
seaty = seat_iy) %>%
select(seatx,seaty,time)
dat_seat## # A tibble: 117 × 3 ## seatx seaty time ## <int> <int> <dbl> ## 1 12 16 2.2 ## 2 12 14 0.8 ## 3 11 5 0.6 ## 4 12 5 2.9 ## 5 5 8 6.2 ## 6 12 15 0.1 ## 7 7 19 14.8 ## 8 9 18 0.7 ## 9 10 18 1.9 ## 10 11 16 0.3 ## # ℹ 107 more rows
◯◯◯ 特定の参加者周辺の近傍同期効果の検証
以下では、特定の座席[6,16]の参加者に関して、近傍同期効果を検証します。 データフレームより、座席[6,16]の参加者のシャッター時間は1.2秒であることがわかります。
dat616 = dat_seat %>% filter(seatx==6 & seaty==16)
dat616## # A tibble: 1 × 3 ## seatx seaty time ## <int> <int> <dbl> ## 1 6 16 1.2
近傍同期効果を具体的に検証するまえに、簡単に近傍同期に関する指標を以下のように定義します。
- [近傍集団] 特定の座席から2席以内の距離にある集団を、その座席にとって「NEIGHBOR」、それ以外の集団を「REMOTE」とする。
- [同期／非同期] 特定の参加者ペアに関して、シャッタータイミングのずれが1秒未満の場合、それらは「SYNC」し、それ以外を「ASYNC」とする。
- [同期率] NEIGHBORまたはREMOTE集団の中で、「SYNC」している参加者の割合を「同期率」と呼ぶ。
以上を踏まえて、座席[6,16]の参加者に関して、NEIGHBORとREMOTEの集団の同期率を比較します。統計的には、NEIGHBORの同期率が有意に高ければ、近傍同期効果は支持されることになります。
# 座席の変数の代入
sx = dat616$seatx;
sy = dat616$seaty;
t = dat616$time
# シャッター時間の差が1秒未満のペアを「SYNC」とし、
# NEIGHBORグループ（周囲2席以内）とREMOTEグループ（それ以外）でSYNCの割合を比較
dat_sync = dat_seat %>%
filter(!(seatx==sx & seaty ==sy)) %>% #同じ席の行（自分自身）を除く
mutate(
# 2席分以内の距離を"NEIGHBOR"、それ以外を"REMOTE"
neighbor = if_else(
(seatx>=sx-2 & seatx<=sx+2) &
(seaty>=sy-2 & seaty<=sy+2),"NEIGHBOR","REMOTE"),
# 1秒未満の距離を"SYNC"、それ以外を"ASYNC"
sync = if_else(
abs(time-t)<1.0, "SYNC", "ASYNC")
)
# 71〜80行目を確認。
dat_sync[71:80,]## # A tibble: 10 × 5 ## seatx seaty time neighbor sync ## <int> <int> <dbl> <chr> <chr> ## 1 5 15 0.3 NEIGHBOR SYNC ## 2 6 14 2.6 NEIGHBOR ASYNC ## 3 4 6 9.8 REMOTE ASYNC ## 4 4 19 9.5 REMOTE ASYNC ## 5 4 5 1.2 REMOTE SYNC ## 6 7 5 10 REMOTE ASYNC ## 7 5 14 10.5 NEIGHBOR ASYNC ## 8 11 12 1.9 REMOTE SYNC ## 9 10 12 9.3 REMOTE ASYNC ## 10 8 10 0.6 REMOTE SYNC
neighborとsyncをグループ化して、2x2の条件ごとに、生起数をカウントします。
dat_sync %>%
group_by(neighbor,sync) %>% #neighbor、syncでグループ分け
count(neighbor, sync, name = "n") #summarizeでn()するのと同様## # A tibble: 4 × 3 ## # Groups: neighbor, sync [4] ## neighbor sync n ## <chr> <chr> <int> ## 1 NEIGHBOR ASYNC 6 ## 2 NEIGHBOR SYNC 7 ## 3 REMOTE ASYNC 67 ## 4 REMOTE SYNC 36
この結果は、座席[6,16]の参加者の近傍グループ（13人）、非近傍グループ（103人）のうち、シャッターのタイミングがSYNCした人は、7人（近傍）および36人 （非近傍）であり、近傍グループの方がシャッターが同調した人が多い（53.8% v.s. 40.0%）ことを示しています。
この結果は、近傍同期効果を支持するものです。
◯◯◯ 全ての参加者に関する近傍同期効果の検証
座席[6,16]に関して検討した解析を、全ての座席について調べてみましょう。 まずは先に座席[6,16]に対して行った処理を、一般の座席に適用させるかたちで関数化します。
compareSyncrate = function(sx,sy,t){
dat_sync = dat_seat %>%
filter(!(seatx==sx & seaty ==sy)) %>% #同じ席の行（自分自身）を除く
mutate(
# ２席分以内の距離を"NEIGHBOR"、それ以外を"REMOTE"
neighbor = if_else(
(seatx>=sx-2 & seatx<=sx+2) &
(seaty>=sy-2 & seaty<=sy+2),"NEIGHBOR","REMOTE"),
# シャッター時間が１秒以内を"SYNC"、それ以外を"ASYNC"
sync = if_else(
abs(time-t)<1.0, "SYNC", "ASYNC")
) %>%
group_by(neighbor,sync) %>%
count(neighbor, sync, name = "n") %>%
ungroup() %>% #以下は、いずれかの条件でn=0のときのための対処
complete( #存在しない組み合わせについても、完全にリストを作成する
neighbor = c("NEIGHBOR","REMOTE"),
sync = c("SYNC","ASYNC"),
fill = list(n = 0) #存在しないcountを0で埋める
)
# いずれかの総数が0の場合、NAを返す。（rateの計算で分母が0になるため）
if((dat_sync$n[1] + dat_sync$n[2] == 0) |
(dat_sync$n[3] + dat_sync$n[4] == 0)){
return(NA)
}
ns = dat_sync$n #総数ベクトル(2x2)
ratio1 = ns[2] / (ns[1] + ns[2]) #NEIGHBOER集団のSYNC率
ratio2 = ns[4] / (ns[3] + ns[4]) #REMOTE集団のSYNC率
return(c(ratio1,ratio2)) #長さ2のベクトルとして返します。
}それでは、実際に座席[6,16]に対する関数の挙動を確認しましょう（正しく計算できていることがわかります）。
dat616 = dat_seat %>% filter(seatx==6 & seaty==16)
compareSyncrate(dat616$seatx,dat616$seaty,dat616$time)## [1] 0.5384615 0.3495146
ひきつづき、
mutateでの使用を想定して、長さ2のベクトルを返すcompareSyncrateを、長さ１の値が返ってくるように関数を少し拡張します。具体的には、新たに引数MODEを立てて、近傍グループの同期率と非近傍グループの同期率を個別に算出できるようにします。
compareSyncrateMode = function(sx,sy,t,MODE){
if(MODE=="NEIGHBOR"){
return(compareSyncrate(sx,sy,t)[1])
}
if(MODE=="REMOTE"){
return(compareSyncrate(sx,sy,t)[2])
}
}
mutateを使って、新たな変数にNEIGHBOR集団とREMOTE集団の同期率を計算していきます。この中で、mapplyを使って、引数をベクトル対応にしている点に注意してください。
dat_seat = dat_seat %>% mutate(
ratioNeighbor = mapply(compareSyncrateMode, seatx,seaty,time,"NEIGHBOR"),
ratioRemote = mapply(compareSyncrateMode,seatx,seaty,time,"REMOTE")
) %>%
filter(!is.na(ratioNeighbor)) #同期率がNAのものを削除
dat_seat## # A tibble: 116 × 5 ## seatx seaty time ratioNeighbor ratioRemote ## <int> <int> <dbl> <dbl> <dbl> ## 1 12 16 2.2 0.25 0.176 ## 2 12 14 0.8 0.333 0.364 ## 3 11 5 0.6 0 0.373 ## 4 12 5 2.9 0.6 0.144 ## 5 5 8 6.2 0.273 0.0571 ## 6 12 15 0.1 0.375 0.259 ## 7 7 19 14.8 0 0.00952 ## 8 9 18 0.7 0.3 0.358 ## 9 10 18 1.9 0.222 0.224 ## 10 11 16 0.3 0.364 0.314 ## # ℹ 106 more rows
グラフにする前に、NEIGHBORとREMOTEのそれぞれの同期率の平均を計算してみましょう。 平均値では、REMOTE集団の19.6%に対して、NEIGHBOR集団は22.5%で、近傍集団では、3%ほど同期率が高いことがわかります。
dat_seat %>%
summarize(
avgNeighbor = mean(ratioNeighbor),
avgRemote = mean(ratioRemote),
)## # A tibble: 1 × 2 ## avgNeighbor avgRemote ## <dbl> <dbl> ## 1 0.225 0.196
（発展）「対応のあるt検定」を行うことで、この差が意味のあるものかどうかを検定することができます。検定の結果、
p<0.05で、有意に、近傍集団の同期率が高いことが示されました。この結果もまた、近傍同期効果を強く支持するものです。
t.test(dat_seat$ratioNeighbor,dat_seat$ratioRemote, paired = TRUE)## ## Paired t-test ## ## data: dat_seat$ratioNeighbor and dat_seat$ratioRemote ## t = 2.109, df = 115, p-value = 0.03712 ## alternative hypothesis: true mean difference is not equal to 0 ## 95 percent confidence interval: ## 0.001791265 0.057152802 ## sample estimates: ## mean difference ## 0.02947203
◯◯◯ Graph6A（NEIGHBOR vs REMOTEの同期率の比較）
# 二つの列変数を行に展開
dat_seat2 = dat_seat %>%
pivot_longer(
c(ratioNeighbor,ratioRemote),
names_to = "group",
values_to = "ratio"
)
# [ggplot]
## x軸をdirection、y軸にdif, 塗りつぶしをdirectionに設定
gp_a = ggplot(dat_seat2,aes(x=group,y=ratio,fill=group))
# [geom_boxplot]ボックスプロットによる描画
gp_a = gp_a + geom_boxplot()
# X軸のラベルの書き換え
gp_a = gp_a + scale_x_discrete(
labels = c("ratioNeighbor"="NEIGHBOR", "ratioRemote"="REMOTE")
)
gp_a = gp_a + scale_y_continuous()
# Bluesパレットを使って、fillを描画。
gp_a = gp_a + scale_fill_brewer(
palette = "Blues",
labels = c("ratioNeighbor"="NEIGHBOR", "ratioRemote"="REMOTE")
)
gp_a = gp_a + labs(x = "TYPE", y = "Ratio",fill="TYPE")gp_a + theme_minimal() +
ggtitle("Graph6A_NeighborSyncEffect")
◯◯◯ Graph6B（NEIGHBOR vs REMOTEの同期率の比較）
Graph6Aは、SYNC率を1秒未満と定義しましたが、これをt秒未満と一般化します。このSYNC率の窓を操作することで、距離による同期効果が、どのような時間スケールで生じるかどうかを検証することができます。
#第４引数に時間窓を追加
compareSyncrateWindow = function(sx,sy,t,t_window){
dat_sync = dat_seat %>%
filter(!(seatx==sx & seaty ==sy)) %>% #同じ席の行（自分自身）を除く
mutate(
# ２席分以内の距離を"NEIGHBOR"、それ以外を"REMOTE"
neighbor = if_else(
(seatx>=sx-2 & seatx<=sx+2) &
(seaty>=sy-2 & seaty<=sy+2),"NEIGHBOR","REMOTE"),
# シャッター時間が１秒以内を"SYNC"、それ以外を"ASYNC"
sync = if_else(
abs(time-t)<t_window, "SYNC", "ASYNC")
) %>%
group_by(neighbor,sync) %>%
count(neighbor, sync, name = "n") %>%
ungroup() %>% #以下は、いずれかの条件でn=0のときのための対処
complete( #存在しない組み合わせについても、完全にリストを作成する
neighbor = c("NEIGHBOR","REMOTE"),
sync = c("SYNC","ASYNC"),
fill = list(n = 0) #存在しないcountを0で埋める
)
# いずれかの総数が0の場合、NAを返す。（rateの計算で分母が0になるため）
if((dat_sync$n[1] + dat_sync$n[2] == 0) |
(dat_sync$n[3] + dat_sync$n[4] == 0)){
return(NA)
}
ns = dat_sync$n #総数ベクトル(2x2)
ratio1 = ns[2] / (ns[1] + ns[2]) #NEIGHBOER集団のSYNC率
ratio2 = ns[4] / (ns[3] + ns[4]) #REMOTE集団のSYNC率
return(c(ratio1,ratio2)) #長さ2のベクトルとして返します。
}
# 長さ１の値を返すように関数を分解
compareSyncrateWindowMode = function(sx,sy,t,t_window,MODE){
if(MODE=="NEIGHBOR"){
return(compareSyncrateWindow(sx,sy,t,t_window)[1])
}
if(MODE=="REMOTE"){
return(compareSyncrateWindow(sx,sy,t,t_window)[2])
}
}
dat_seatを初期化します。
dat_seat = dat.fs %>%
filter(exp=="Single") %>%
filter(!is.na(seat)) %>% #seatが不明なものを取り除く。
mutate(seatx = seat_ix, #変数名を変えます。
seaty = seat_iy) %>%
select(seatx,seaty,time)
dat_seat## # A tibble: 117 × 3 ## seatx seaty time ## <int> <int> <dbl> ## 1 12 16 2.2 ## 2 12 14 0.8 ## 3 11 5 0.6 ## 4 12 5 2.9 ## 5 5 8 6.2 ## 6 12 15 0.1 ## 7 7 19 14.8 ## 8 9 18 0.7 ## 9 10 18 1.9 ## 10 11 16 0.3 ## # ℹ 107 more rows
getSyncWindow = function(t_window){
dat_seat %>%
mutate(
ratioNeighbor = mapply(
compareSyncrateWindowMode, seatx,seaty,time,t_window,"NEIGHBOR"),
ratioRemote = mapply(
compareSyncrateWindowMode,seatx,seaty,time,t_window,"REMOTE")
) %>%
filter(!is.na(ratioNeighbor)) %>%
summarize(
avgNeighbor = mean(ratioNeighbor),
avgRemote = mean(ratioRemote)
)
}
getSyncWindowMode = function(t_window,MODE){
avgs = getSyncWindow(t_window)
if(MODE=="NEIGHBOR"){
return(avgs$avgNeighbor)
}else{
return(avgs$avgRemote)
}
}
# t_window=1.0のときの結果を確認（正しいですね）
getSyncWindowMode(1.0,"NEIGHBOR"); getSyncWindowMode(1.0,"REMOTE")## [1] 0.2253262
## [1] 0.1958542
#mapply(getSyncWindow,0.2 * (1:25))
dat_wi = tibble(
#t_window = 0.1 * (0:15)
t_window = 0.2 * (0:10)
#t_window = 1.0
)
dat_wi = dat_wi %>%
mutate(ratioNeighbor =
mapply(getSyncWindowMode,t_window,"NEIGHBOR"),
ratioRemote =
mapply(getSyncWindowMode,t_window,"REMOTE")
) %>%
pivot_longer(
cols = c(ratioNeighbor,ratioRemote),
names_to = "group",
values_to = "ratio"
)
dat_wi## # A tibble: 22 × 3 ## t_window group ratio ## <dbl> <chr> <dbl> ## 1 0 ratioNeighbor 0 ## 2 0 ratioRemote 0 ## 3 0.2 ratioNeighbor 0.0491 ## 4 0.2 ratioRemote 0.0562 ## 5 0.4 ratioNeighbor 0.102 ## 6 0.4 ratioRemote 0.103 ## 7 0.6 ratioNeighbor 0.160 ## 8 0.6 ratioRemote 0.143 ## 9 0.8 ratioNeighbor 0.204 ## 10 0.8 ratioRemote 0.174 ## # ℹ 12 more rows
gp_b = ggplot(dat_wi,aes(x=t_window,y=ratio,colour=group))
gp_b = gp_b + geom_line(size=0.5) + geom_point(size=3)
gp_b = gp_b + scale_colour_brewer(
palette ="Set2",
labels = c("ratioNeighbor"="NEIGHBOR", "ratioRemote"="REMOTE")
)
gp_b = gp_b + scale_x_continuous(limits=c(0,2.0))
gp_b = gp_b + scale_y_continuous(limits=c(0,0.4))
gp_b = gp_b + labs(x="TIME WINDOW[s]",y="SYNC RATE",colour="GROUP TYPE")gp_b + theme_minimal() +
ggtitle("Graph6B_NeighborSyncEffect")
グラフから、近傍同期効果が現れるのは、およそ0.5秒以降であることがわかります。ちなみに、実際に検定をしてみるとNEIGHBORとREMOTEとの間に有意な差があるのは、TIME WINDOWが1.0秒付近の限られた区間のみとなります。つまり、距離的に近いもの同士が影響を与え合う時間窓は、およそ1秒程度であることが示唆されます。
[確認課題]
ApMedia06_Work
- 20秒間で10回押すフリーシャッター実験の結果を用いて、「コンマN秒」の分布を作成していください（使用するデータフレームは以下の解説を参照）。
- この際、前半1〜5回目のシャッターを「1ST」、後半6〜10回目のシャッターを「2ND」として、前後半別で分布を重ねて示してください。
- 締切は11月中とします. -ファイル名は「2350xx_work6.R」としてください。
[10times フリーシャッター実験]
20秒の間にシャッターを10回押してください実験の結果（24年7月11日）が以下から入手できます。
url = "https://lab.kenrikodaka.com/_download/csv/freeshutter_240711.csv"
dat_240711 = read.csv(url)参加者は109人
dim(dat_240711)## [1] 109 12
最初と最後の6行
dat_240711[c(1:6,104:109),]## seatx seaty t1 t2 t3 t4 t5 t6 t7 t8 t9 t10 ## 1 9 2 0.3 1.1 1.4 1.6 1.7 2.1 2.6 2.7 3.6 5.1 ## 2 12 9 0.4 1.2 1.6 2.3 3.0 13.3 14.0 15.0 16.2 19.0 ## 3 4 4 1.2 2.9 4.9 6.8 9.0 10.7 12.5 15.0 17.2 19.2 ## 4 10 2 6.0 6.5 6.6 6.9 7.2 7.5 7.7 7.9 8.2 8.4 ## 5 6 14 1.3 2.7 3.2 8.2 8.9 11.2 13.3 17.3 17.7 18.0 ## 6 9 6 1.1 1.5 2.1 2.4 3.5 4.0 5.4 7.7 9.8 11.1 ## 104 9 1 1.1 2.4 3.8 5.5 6.9 8.3 9.2 10.8 12.1 NA ## 105 12 4 1.2 2.4 4.2 5.2 6.2 7.4 8.2 9.4 10.5 11.5 ## 106 12 15 2.3 3.4 4.5 6.7 10.3 11.8 16.9 18.9 20.1 21.0 ## 107 10 15 0.9 1.6 2.2 3.5 4.6 5.2 5.9 6.5 8.1 9.9 ## 108 12 15 2.0 3.6 4.6 5.7 10.2 11.5 16.7 17.9 18.7 19.9 ## 109 11 2 2.9 3.1 3.6 4.0 5.2 14.2 15.7 16.0 16.7 18.8
- seatx：座席の水平方向（生徒視点で左から1,2,…..）
- seaty：座席の奥行き方向（教卓に近い方から1,2,…)
- t1, t2, t3, …：1回目のシャッター時間, 2回目のシャッター時間, …
- NAは欠損データ（例えば、104行目のt10）