[R]iPhoneのヘルスケアの歩数データを分析

平成最後に、iPhoneのアプリ「ヘルスケア」のデータをダウンロードして、自分の歩数データで遊んでみました。

#ヘルスケアアプリデータのダウンロード
- [参考：iPhoneに記録された「ヘルスケア」データを書き出すには?](https://news.mynavi.jp/article/20170910-ipadiphonehacks/)

>iPhoneには、万歩計相当の機能が内蔵されている。iPhone 5s以降の機種には、運動に関する情報をわずかな消費電力で測定できるモーション・コプロセッサが標準装備されているため、特になにもしなくても歩数や移動距離を記録できるのだ。


かれこれ4-5年はiPhoneユーザでしたが、データを確認してみるといい感じで蓄積されていたので、集計して遊んでみました。
自分の歩数、移動距離にどんな傾向があるのかと考えてみましたが、

- 平日は仕事で基本は会社と家の往復。天気が悪くても同じような傾向かも。
- 土日祝日は外出して極端に遠出するか、家を出ない・近場のみの移動。雨の日はあんまり家を出ないかも。
- 昔より運動しなくなったから歩幅とかちょっと狭まっているかも。
- 引っ越し後、会社が近くなったので、移動距離が減っているかも。

正しいかどうかは置いておいて、蓄積されたデータから見て取れるのか試してみます。

#データのcsv化
ダウンロードしたデータは以下のようなxmlファイルとなっています。
適当なブラウザやメモ帳等で確認。
Recordタグの中に、身長・体重、歩数や日付のデータが格納されていることが分かります。
分析に使えそうな以下のデータを取り出したいです。

1. value : 歩数または移動距離の数値データ
2. type : 歩(~StepCount)または移動距離(~DistanceWalkingRunning)のカテゴリデータ
3. ~Date : 時間の区切りのタイミングは不明ですが、カウント時の日付＋時間データでしょう。

<img src="https://i.imgur.com/yGDhtkS.png" width="100%">


csvファイル等のテーブルデータへの変換の方法はいろいろあるようですが、ここではRのXMLライブラリを使って変換してみます。


- [参考：R による XML の CSV 化](http://fits.hatenablog.com/entry/2013/09/23/002707)



```r:csvへの変換
library(XML)
library(lubridate)

doc <- xmlParse("healthcare.xml")
items <- getNodeSet(doc, "//Record")

#value/startDate/typeの書き出し
value <- sapply(items, function(x) xmlGetAttr(x, "value"))
startdate <- sapply(items, function(x) xmlGetAttr(x, "startDate"))
type <- sapply(items, function(x) xmlGetAttr(x, "type"))

#データフレーム化
df <- data.frame(value, startdate, type, stringsAsFactors = FALSE)

#データ確認
> head(df,10)
   value                 startdate                              type
1    *** 2017-03-08 06:29:20 +0900    HKQuantityTypeIdentifierHeight
2     ** 2017-03-08 06:29:20 +0900  HKQuantityTypeIdentifierBodyMass
3     46 2018-11-27 09:59:21 +0900 HKQuantityTypeIdentifierStepCount
4      3 2018-11-27 10:16:59 +0900 HKQuantityTypeIdentifierStepCount
5      8 2018-11-27 11:47:56 +0900 HKQuantityTypeIdentifierStepCount
6    128 2019-02-10 15:39:55 +0900 HKQuantityTypeIdentifierStepCount
7     70 2019-02-10 15:51:48 +0900 HKQuantityTypeIdentifierStepCount
8     26 2019-02-10 16:17:48 +0900 HKQuantityTypeIdentifierStepCount
9     13 2019-02-10 16:26:56 +0900 HKQuantityTypeIdentifierStepCount
10    23 2015-03-05 06:24:08 +0900 HKQuantityTypeIdentifierStepCount

```

#分析用データの作成
取り出した3列のデータを使って、分析に使いたいテーブルを作成してみます。
冒頭に書いた予想に必要になりそうな情報を取り入れてみます。

①平日と休日の違いを見たいので、祝日データを拝借して追加。 [参考：祝日のデータ]("https://github.com/takitsuba/yeah_meccha")
②個人的な自分の休みの日（有休とか）を追加
③気象庁から東京の降水量のデータを拝借して追加。[参考：過去の気象データ検索]("https://www.data.jma.go.jp/obd/stats/etrn/index.php")

```r
library(magrittr)
library(tidyverse)
library(ggplot2)


#日付データのみ取り出し日付型に、valueを数値型に
df$date <- unlist(str_extract_all(df$startdate, "\\d{4}-\\d{2}-\\d{2}")) %>% as.Date(.)
df$value <- as.numeric(df$value)


#歩数データと移動距離データの作成
df_step <- df %>% filter(type == "HKQuantityTypeIdentifierStepCount")
df_km <- df %>% filter(type == "HKQuantityTypeIdentifierDistanceWalkingRunning")

#日付ごとに歩数・移動距離を集約
df_step2 <- df_step %>% group_by(date) %>% summarise(stepcount = sum(value))
df_km2 <- df_km %>% group_by(date) %>% summarise(km = sum(value))

#連続日付データを作成
date <- seq(as.Date("2015-03-12"), as.Date("2019-04-11"), by = "day")
date <- as.data.frame(date)
colnames(date) <- "date"

#日付をキーに、歩数と移動距離データを結合
df_1 <- date %>% left_join(df_step2, by = "date") %>% 
  left_join(df_km2, by ="date")

#曜日データ作成
df_1$wday <- wday(df_1$date, label = TRUE)
#歩幅のデータを作成
df_2 <- df_1 %>% mutate(length_m = km/stepcount*1000)


#祝日の表示
holiday <- read.csv("holiday.csv", header = F)
colnames(holiday) <- c("date", "holiday")
holiday$date <- as.Date(holiday$date)
holiday$holiday <- str_detect(holiday$holiday, ".")


#有休の反映（個人的にカレンダーを追って作成）
add_holiday <- read.csv("kyuka.csv", header = T)
add_holiday$date <- as.Date(add_holiday$date)

#降水量の反映
rain <- read.csv("rain.csv", header = T)
rain$date <- as.Date(rain$date)

#各種外部データを日付をキーに結合
df_3 <- df_2 %>% left_join(holiday, by = "date") %>% 
  left_join(add_holiday, by = "date") %>% 
  left_join(rain, by = "date")  
  

#平日と休みをweekdayとholidayの2パターンに変換
df_3$wday2[df_3$wday %in% c("月","火","水","木","金")] <- "weekday"
df_3$wday2[df_3$wday %in% c("土","日")] <- "holiday"
df_3$wday2[df_3$holiday == TRUE] <- "holiday"
df_3$wday2[df_3$add_holiday == TRUE] <- "holiday"
df_3$wday2 <- as.factor(df_3$wday2)
 
df_3 %<>% select(-c(holiday, add_holiday))
```
上のしょっぱい操作で以下のようなcolumnを持つテーブルを作成しました。
「date:日付」「stepcount:歩数」「km:移動距離(km)」「wday:曜日」「length_m:歩幅(m)」「rain:降水量(mm)」「wday2:平日か休日の分類」
データを確認してみます。

```r
#データの確認
> head(df_3)
        date stepcount      km wday  length_m rain   wday2
1 2015-03-12      7809 5.42587   木 0.6948226  0.0 weekday
2 2015-03-13      9395 6.65875   金 0.7087547  0.0 weekday
3 2015-03-14      7965 5.20370   土 0.6533208  0.0 holiday
4 2015-03-15      5610 4.01634   日 0.7159251  0.0 holiday
5 2015-03-16     10198 7.16307   月 0.7023995 10.0 weekday
6 2015-03-17      8891 6.47283   火 0.7280205  0.5 weekday


#NAの有無の確認
> sapply(as.data.frame(df_3), function(x) any(is.na(x)))
     date stepcount        km      wday  length_m      rain     wday2 
    FALSE      TRUE      TRUE     FALSE      TRUE     FALSE     FALSE 

#サマリーの確認
> summary(df_3)
      date              stepcount           km           wday        length_m           rain             wday2    
 Min.   :2015-03-12   Min.   :   34   Min.   : 0.03584   日:213   Min.   :0.3673   Min.   :  0.000   holiday:551  
 1st Qu.:2016-03-18   1st Qu.: 6008   1st Qu.: 4.05888   月:213   1st Qu.:0.6434   1st Qu.:  0.000   weekday:941  
 Median :2017-03-26   Median : 7172   Median : 4.83882   火:213   Median :0.6781   Median :  0.000                
 Mean   :2017-03-26   Mean   : 7368   Mean   : 4.97775   水:213   Mean   :0.6729   Mean   :  4.307                
 3rd Qu.:2018-04-03   3rd Qu.: 8758   3rd Qu.: 5.99860   木:214   3rd Qu.:0.7067   3rd Qu.:  1.500                
 Max.   :2019-04-11   Max.   :26988   Max.   :18.23305   金:213   Max.   :1.4503   Max.   :156.500                
                      NA's   :43      NA's   :43         土:213   NA's   :43  

```

43日分欠損データがあるようですが、サマリーから歩数他ざっくり、概要が分かりました。
一日平均
歩数：7400歩
移動：5km
歩幅：67cm

平日は家と会社の往復のみだけど、意外と歩いてるなーという印象。
歩幅に関しては成人男性の平均は70cmくらいらしいので、若干の短足？が疑われます。



#データの確認
作成したテーブルの概要をもう少し把握するため、データを可視化してみます。
ペアプロット作成のために`GGally`の`ggpair`関数を使用しました。

- [参考：GGallyパッケージのggpair関数を使いこなすための覚え書き](https://assabu.exblog.jp/28933942/)

```r
library(GGally)
#pairplotの作成
ggpairs(df_3, aes(color = wday2))
```

<img src="https://i.imgur.com/v48UIOp.png" width="100%">

データ概要を見てみると、思った通り平日より休日のほうが歩数や移動距離のばらつきが大きそうです。
一方、降水量はこのままだとあまり影響のほどが読み取れないので、降水量データを「雨」か「晴れ」かの2カテゴリに変換してみます。
一般的に、降水量3mm以上だと「普通に雨が降っている状況」らしいので、降水量の列で3mmより大きい数値を「rain」、それ以下を「sunny」に変換した列を作成してみます。

- [参考：降水量の目安のイメージ！1mmから100mmまでを調べてみた](https://happy-teacher-singing.com/archives/2651)

```r
#3mmを基準に雨のデータをカテゴライズした列を作成
df_3 %<>% mutate(rain2 = as.factor(ifelse(rain > 3, "rain", "sunny")))

#データの確認
> head(df_3)
        date stepcount      km wday  length_m rain   wday2 rain2
1 2015-03-12      7809 5.42587   木 0.6948226  0.0 weekday sunny
2 2015-03-13      9395 6.65875   金 0.7087547  0.0 weekday sunny
3 2015-03-14      7965 5.20370   土 0.6533208  0.0 holiday sunny
4 2015-03-15      5610 4.01634   日 0.7159251  0.0 holiday sunny
5 2015-03-16     10198 7.16307   月 0.7023995 10.0 weekday  rain
6 2015-03-17      8891 6.47283   火 0.7280205  0.5 weekday sunny
```
では実際、「休日と平日」、「雨が降っているか降っていないか」で歩数の平均に差があるのか、見てみます。
t検定を休日グループと平日グループのそれぞれで「雨」「晴れ」の場合分けで実施してみました。

```r
#「休日」の雨、晴れで検定
df_3_holiday_rain <- df_3 %>% filter(wday2 == "holiday" & rain2 == "rain") %>% select(stepcount)
df_3_holiday_sunny <- df_3 %>% filter(wday2 == "holiday" & rain2 == "sunny") %>% select(stepcount)

> t.test(df_3_holiday_rain, df_3_holiday_sunny)
	Welch Two Sample t-test

data:  df_3_holiday_rain and df_3_holiday_sunny
t = -2.0834, df = 142.33, p-value = 0.039
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -2027.26553   -53.23771
sample estimates:
mean of x mean of y 
 6097.208  7137.460 


#「平日」の雨、晴れで検定
df_3_weekday_rain <- df_3 %>% filter(wday2 == "weekday" & rain2 == "rain") %>% select(stepcount)
df_3_weekday_sunny <- df_3 %>% filter(wday2 == "weekday" & rain2 == "sunny") %>% select(stepcount)

> t.test(df_3_weekday_rain, df_3_weekday_sunny)

	Welch Two Sample t-test

data:  df_3_weekday_rain and df_3_weekday_sunny
t = 0.15803, df = 336.28, p-value = 0.8745
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -277.3411  325.7957
sample estimates:
mean of x mean of y 
 7630.146  7605.918 
```

結果を見てみると、平日の場合はp値0.8745、休日の場合はp値0.039ということで、
平日は雨か晴れかは歩数に与える影響はあまりなく、休日の場合は大いにあると言えそう。

最初の予想がデータからも読み取れていると言えそうです。
改めて、場合分けで可視化してみると、やはり違いがあり、休日の雨の日は、やっぱり外出は控える傾向にあると言えます。

```R
ggplot(df_3, aes(x = rain2, y = stepcount, fill = rain2)) +
  geom_boxplot()  + 
  theme(legend.position = 'none') +
  facet_wrap(~wday2)
```
<img src="https://i.imgur.com/UKr11dB.png" width="100%">


時系列の歩数変化も少し見てみます。

```r
ggplot(df_3, aes(x = date, y = stepcount, colour = wday2)) + 
  geom_point() +
  geom_smooth()
```


<img src="https://i.imgur.com/LXCd4GR.png" width="100%">


#モデルの作成
最後に、これまでの結果から、「休日か平日か」「雨か晴れか」で歩数に影響がありそうなことが分かったので、その日の天候と自分の歩数から、「休日か平日か」を予測する決定木モデルを作ってみました。
（本質的にあまり意味のない分析であることはご了承ください。）

```r:決定木による分類
library(rpart)
library(rattle)
library(rpart.plot)

#2018年までをトレーニングデータ
df_train <- df_3 %>% filter(date >= "2015-03-12" & date <= "2018-12-31") %>% select(wday2, stepcount, rain2)
#2019年からを予測用データに分割
df_pred <- df_3 %>% filter(date >= "2019-01-01" & date <= "2019-04-11") %>% select(wday2, stepcount, rain2)

#テストデータをもとにモデル作成。cpで分類の細かさを指定
tree <- rpart(wday2~stepcount + rain2, df_train, cp = 0.012)
#モデルの描画
fancyRpartPlot(tree)
```


<img src="https://i.imgur.com/xNdypHd.png" width="100%">


```r:モデルの検証
pred <- predict(tree, df_pred)

#予測用データに結果を結合し、正解率算出
result <- cbind(df_pred, pred)
result$pred <- ifelse(result$holiday >= result$weekday, "holiday", "weekday")
result$answer <- ifelse(result$wday == result$pred, 1, 0)
answer <- sum(result$answer)/length(result$answer)

> answer
[1] 0.7623762
```
上記のモデルで76%程度の予測ができました。

#おわりに
自分の歩数データを使用して遊んでみました。
天気や休日のデータ以外にももっと面白い考察もできそうです。
今後もいろんなデータでコツコツ、粛々と理解を深めていきたいところです。