歩数データから行動パターンを抽出する話

この1年ほどfitbitで歩数データを取り続けている。

最近、歩数データは1日よりも小さい単位（1分、15分）で取得できることを知ったので、それを利用して日々の行動をパターン化したい。

同時に手持ちの睡眠データと組み合わせて、行動パターンとの関連も見たい。

なお、2016年に入ってからジムに行き始めたりいろいろと行動が変わったこともあるので、対象期間は本日4/21までの90日間とする。

パッケージの読み込み

library("myFitbit")
library("dplyr")
library("lubridate")
library("httr")
library("jsonlite")
library("topicmodels")
library("tidyr")
library("slam")
library("openair")
library("ggplot2")
library("gridExtra")
library("readxl")

fitbitからデータを取得する。

fitbit APIを用いてデータを取得する。

fitbitrパッケージという便利なパッケージもあるが、いろいろとカスタマイズしたかったので自分のパッケージmyFitbitを用いた。

一日の定義を自身の生活時間に合わせて午前4時から24時間としていることに注意。

getStepHourlydata <- function(date, token){
  date2 <- as.Date(date) + 1
  u <- sprintf("https://api.fitbit.com/1/user/-/activities/steps/date/%s/%s/15min/time/04:00/03:00.json",
               date, date2)
  res <- httr::GET(url = u,
                   config(token = token))
  dat <- fromJSON(content(res, as = "text"))
  result <- dat$`activities-steps-intraday`$dataset
  result$time <- as.POSIXct(paste(date, result$time))
  result <- result %>% group_by(hour=ifelse(hour(time) >= 4, hour(time), hour(time) + 24)) %>% summarise(steps = sum(value))
  result$date <- date
  return(as.data.frame(result, stringsAsFactors = FALSE))
}

fitbit_token <- getToken("MY KEY", "MY SECRET")
result <- bind_rows(lapply(seq(Sys.Date()-90, Sys.Date(), by = 1), getStepHourlydata, token = fitbit_token))
result_wide <- result %>% mutate(hour = paste0("H", formatC(hour, width = 2, flag = "0"))) %>% spread(hour, steps)

LDAでクラスタリングする。

この研究会でLDAを歩数パターンに適用している発表があったので、追試の意味をこめてLDAでクラスタリングする。

平日、休日で2パターンずつ、あとその他グループが1つくらい見つかればいいかなという考えで、グループ数は5とする。

LDAにはtopicmodelsパッケージのLDA関数を用いる。

この際、入力データがdata.frameだと受け付けてくれないので、slamパッケージの関数を用いて、simple_triplet_matrixに変換しておく。

# LDA modeling
result_stm <- as.simple_triplet_matrix(as.data.frame(result_wide[,-1]))
model_lda <- LDA(result_stm, 5)
inf_lda <- posterior(model_lda, result_stm)

tbl_topics <- table(topics(model_lda))

result_lda_group <- data.frame(inf_lda$terms) %>% mutate(group = paste0("Group ",1:5, " (n = ", tbl_topics, ")")) %>% gather(var, value, -group) %>% mutate(hour = extract_numeric(var))
th_value <- 0.001
result_lda_group$value <- ifelse(result_lda_group$value < th_value, 0, result_lda_group$value)

日単位で歩数を集計し、ローカルに保有している自身の睡眠データと結合する。

# daily data
result_daily <- result %>% group_by(date) %>% summarise(steps = sum(steps))
result_daily$group <- factor(topics(model_lda), levels = 5:1)


# join my sleep data
d_sleep <- read_excel("result_all.xlsx")
d_sleep$date <- as.Date(d_sleep$date)
result_sleep <- d_sleep %>% inner_join(result_daily)

結果を可視化する。

各グループの特徴を時間別の分布、総歩数、睡眠の質（入眠後3時間以内の深い眠りの割合）で確認する。

# Group distribution
p_line <- ggplot(result_lda_group, aes(x = hour, y = value, color = group)) + geom_line() + geom_point() +
  scale_color_brewer(palette = "Set1") +
  scale_x_continuous(breaks = 4:27)+
  facet_wrap(~group, ncol= 1) + theme(legend.position = "none")
# Total steps boxplot
p_bx <- ggplot(result_sleep, aes(x = group, y = steps)) + geom_boxplot(fill = "#edc58f") + coord_flip()
# Sleep boxplot
p_bx_sleep <- ggplot(result_sleep, aes(x = group, y = per_D_180)) + geom_boxplot(fill = "#b4ebfa") + coord_flip()

grid.arrange(
  p_line, p_bx, p_bx_sleep,
  widths = c(2, 1, 1)
)

分類した5つのグループおよび総歩数がどのように分布しているかカレンダープロットで確認する。

こちらはグループをカレンダープロットで可視化したもの。

# group of each days
lattice::trellis.par.set(grid.pars = list(fontfamily = "HiraKakuProN-W3"))
calendarPlot(result_daily, pollutant = "group", year = "2016", cols="Set1", 
             labels = paste0("Group ", 1:5), breaks = c(0:5), w.shift = 2)

こちらは総歩数をカレンダープロットで可視化したもの。

# total steps of each days
calendarPlot(result_daily, pollutant = "steps", year = "2016", w.shift = 2)

知見

• LDAの結果をカレンダープロットで確認すると、Group1は平日パターン、Group3,5は休日パターンのように見える。Group2,4は飲み会パターン、夜間の歩数が多いので飲み会パターンか（弊社はホワイトなので深夜に帰宅することはまず無い）。
• 総歩数はどのパターンも中央値はあまり変わらない。休日パターン（Group5）は左に裾をひいてるので概して歩数が少なくなっている印象。
• 平日パターン（Group1）と休日パターン（Group5）の眠りの質（per_D_180）を比べると、前者の方が良さそう。とはいえ、ばらつきも大きい。総歩数のデータと併せて考えると、良い睡眠に一定の運動量は必要だが、それが十分条件ではなさそう。

知見をふまえて決定木でモデリングする。

上記知見をふまえて決定木でモデリングしてみる。

先の睡眠データには「就寝前2時間以内に食事をしたか(shokuji)」「就寝前1時間以内にブルーライトを発するような装置を注視したか(BL)」「寝る1時間前に蛍光灯をつけたか(Keikou)」がそれぞれ0/1（いいえ/はい）で含まれている。

また、総歩数はstepsとして格納されている。

以上をrpartを用いてモデリングする。

library("partykit")
library("rpart")
result_sleep$shokuji <- factor(result_sleep$shokuji)
result_sleep$BL <- factor(result_sleep$BL)
result_sleep$Keikou <- factor(result_sleep$Keikou)
model_dt <- rpart(data = result_sleep, per_D_180 ~ group + steps + shokuji + BL + Keikou)
plot(partykit::as.party(model_dt))

結果をみるとブルーライトを注視しないと圧倒的に睡眠の質が良い。

またばらつきは大きいものの、歩数パターンのGroup1と3は他のグループに比べて睡眠の質が良い。

Group1は主に平日に現れる歩数パターンなので、平日の方が休日より睡眠の質が良いのか。

また悪い歩数パターンの中でも約8,000歩を超えると多少は睡眠の質が良くなっている。 当面の目標値は8000歩にしておこう。