The 70th Tokyo.R

画像の出展(https://poseidonhd.com/tvshows/violet-evergarden/)

Introduction

本稿は、第70回 Tokyo.R 『ベイズ・統計モデリング』 応用セッションへのエントリーです。応用セッションで紹介させていただいたベイズ統計モデリングの実践例を解説します。発表スライドは非公開ですが、本記事の方が情報量が多いため、その点どうかご容赦くださいませ。

本稿で利用した解析スクリプト、stanファイル、一時保存データなどは、全てGitHubに置いてあります。VioretEvargardenというフォルダです。再現実験やトレーニングにご活用くださいませ。

Topics

  1. ヴァイオレット・エヴァーガーデン(Violet Evergarden)
  2. Twitterフォロワー数の取得(Collect Twitter Data)
  3. ローデータの可視化(Data Visualization)
  4. データ生成メカニズムと統計モデリング(Statistical Modeling)
  5. Stan用のデータ整形(Data Preparation for Stan)
  6. Stanによる統計モデルの実装(Stan Model)
  7. 結果の確認と可視化(Result)
  8. 総括(Conclusion)

Libraries

本稿で使用する各種ライブラリです。事前のインストールをお願い致します。

# libraries
library(RCurl) # API叩く時に
library(rlist) # list操作
library(rjson) # json操作
library(tidyverse) # tidyパッケージたち(主にデータ整形dplyr&tidyrと可視化ggplot2)
library(lubridate) # 時間変数の操作
library(anytime) # 時間変数の操作
library(formattable) # HTML表
library(DT) # HTML表
library(scales) # dat_breaks()でggplot2の時間ラベル
library(shinystan) # stan結果確認と収束診断
library(ggmcmc) # stan可視化とサンプル取り出し
library(rstan) # stan
# stan並列化とコンパイル上書き
options(mc.cores = parallel::detectCores())
rstan_options(auto_write = TRUE)

1. ヴァイオレット・エヴァーガーデン(Violet Evergarden)


本稿では、2018年冬期放送のアニメーション作品『ヴァイオレット・エヴァーガーデン』のTwitterフォロワー数をテーマに取り上げます。本作品は全13話から構成されており、放送期間は2018年の1月10日から同年4月4日でした。京都アニメーション大賞を受賞。絵がとっても綺麗で、作中の人物に心が惹かれます。当期のNo.1アニメだったと言えるかもしれません。個人的には本作品が一番のお気に入りでした。

あらすじや詳細は、ヴァイオレット・エヴァーガーデン公式サイトを御覧ください。


2. Twitterフォロワー数の取得(Collect Twitter Data)


ヴァイオレット・エヴァーガーデンには、公式のTwitterアカウントがあります。フォロワー数7万6065人(2018年5月7日時点)!!スゴイですね。どのようにフォロワー数というのは増えていくのでしょうか。確認したいと思います。

アニメ公式Twitterアカウントのフォロワー数のログを取っている有り難いAPIがありますので、こちらからデータを取得したいと思います。また、本APIをRから操作する際の方法は、いずにゃんの研究日記にまとめられています。あわせてご参照くださいませ。

それでは、ヴァイオレット・エヴァーガーデン公式アカウントフォロワー数のデータ取得から行っていきます。

Master情報の取得

まずは、Master情報を取得します。どんなデータがあって、どのようなkeyで引き出さなければならないのか。確認します。

# Referrence
  # http://izunyan.hatenablog.com/
  # https://github.com/Project-ShangriLa/sora-playframework-scala

#アニメの放送年:ベクトル(複数指定可)
yearsBC<-c(2018)
#アニメシーズン:ベクトル(1=冬、2=春、3=夏、4=秋: 複数指定可能)
seasons<-c(1)
#現在から何時間前までのデータを取得するか
endTime<-c(1)
#50時間 ✕ 何回分のデータを取得するか
samplings<-c(56)
  # as.POSIXct("2018-05-06") - as.POSIXct("2018-01-08") # データ取得日(5月6日)から、APIデータ記録開始の1月8日まで
  # 118 * 24 / 50 # 50時間を1単位で何回取ればよいか

# Masterデータを取得し、Twitter公式アカウント名とアニメ名を取り出す
master<-data.frame(matrix(NA,0,4))
for(i in 1: length(yearsBC)){
  for(j in 1: length(seasons)){
    tmpFull <- paste("http://api.moemoe.tokyo/anime/v1/master",yearsBC[i],seasons[j],sep = "/") %>% 
      sprintf("Array") %>%
      list.load() %>% list.ungroup() 
    N_anime<-length(tmpFull)/15
    for (k in 1:N_anime){
      tmpTwAcName<-unlist(tmpFull[c((k*15)-13)][1])
      tmpAnime<-unlist(tmpFull[c((k*15)-11)][1])
      master[(length(master[,1])+1),]<-c(yearsBC[i],j,tmpTwAcName,tmpAnime)
    }
  }
}

# マスター情報の確認
DT::datatable(master[33:35,])

ヴァイオレット・エヴァーガーデンのアカウント名は33番目あたりでしょうか。“Violet_Letter”というアカウント名を使ってデータを引き出します。

Followerデータの取得

# endpointの決定
subtract_UNIXtime<-c()
for(t in 1:samplings){subtract_UNIXtime[t]<-c((((60*30*100)*t)-(60*30*100))-3)}
end<-c(as.numeric(Sys.time())-c(endTime*60))
endP<-round(c(end-subtract_UNIXtime)-180000,0)

# Follower Historyの抽出(ヴァイオレット・エヴァーガーデン)
FollowerHistory<-data.frame(matrix(NA,0,5))
for(t in 1: samplings){
  tmpRes <- paste("http://api.moemoe.tokyo/anime/v1/twitter/follower/history?account=Violet_Letter&end_date=",endP[t],sep="") %>%
    sprintf("account") %>%
    list.load("") %>% list.ungroup()
  for(u in 1 : (length(tmpRes) / 2 )){
    FollowerHistory[(length(FollowerHistory$X1)+1),]<-c(2018,1,"ヴァイオレットエヴァーガーデン",as.numeric(unlist(tmpRes[[2*u]][1])),as.numeric(unlist(tmpRes[[(2*u)-1]][1])))
  }
}

これで、データの取得自体は完了です。

ローデータの整形

早速フォロワー数の推移を見たいところですが、データの整形が必要です。欠損やUNIX時間が含まれていたりするので、少し回りくどい整え方になっています。放送日の間歇データをjoinするのも、ちょっと手間ですね。読みにくくてすみません。

#データの列名や時間情報を整理
base <- data.frame(Date = as.Date(seq(as.POSIXct("2018-01-09"), as.POSIXct("2018-05-07"), by = "days")))
TidyAnime<- FollowerHistory %>% 
  dplyr::mutate(Year = as.integer(X1),
                Season = as.integer(X2),
                Anime = as.factor(X3),
                UnixTime = as.integer(X4),
                Follower = as.integer(X5)) %>%
  dplyr::mutate(Time = as.POSIXct(as.Date(anytime(as.numeric(UnixTime),tz = "Asia/Tokyo"))),
                Date = lubridate::date(Time)
  ) %>%
  dplyr::select(Year, Season, Anime, Time, Date, Follower) %>%
  dplyr::group_by(Date) %>%
  dplyr::arrange(desc(Time), .by_group = TRUE) %>%
  dplyr::distinct(Date, .keep_all = TRUE) %>%
  dplyr::ungroup() %>%
  dplyr::right_join(base) %>%
  dplyr::full_join(data.frame(Date = lubridate::date(as.POSIXct(c("2018-01-09"))) + (seq(1:13) * 7 - 6),
                              OnAir = rep(1),
                              OnAirDate = lubridate::date(as.POSIXct(c("2018-01-09"))) + (seq(1:13) * 7 - 6))
  ) %>%
  dplyr::arrange(Date) %>%
  tidyr::replace_na(list(OnAir = 0)) %>%
  tidyr::fill(Year, Anime, Season, .direction = "down") %>%
  dplyr::filter(Date > as.POSIXct("2018-01-08") & Date < as.POSIXct("2018-04-18"))%>%
  dplyr::mutate(OnAir = ifelse(OnAir == 1, "放送日", "非放送日")) %>%
  as.data.frame()

# .Rdataで保存
save(TidyAnime, file = "data/TidyAnime.Rdata")

本APIは、指定した期間で取得できるデータの量に制限があります。今回は「今日から50時間前まで」「50時間まえから100時間前まで」…と、データ取得日に依存してしまう取得の方法を採用しています。よって、解析結果を再現できるように取得したデータを.Rdataで保存し、読み込んで使います。

再現性のために.Rdataを共有します。

.Rdataだけを必要とする場合は、下記からダウンロードし、読み込んでください。

TidyAnime.Rdataのダウンロード先(うなどんGitHub)

データの確認

保存したデータを読み込みます。

画面の都合上、いくつかの変数に絞って提示します。これで、ベースのデータセットが整理出来ました。

load(file = "data/TidyAnime.Rdata")
DT::datatable(TidyAnime %>% dplyr::select(1:3, Date, Follower, OnAir))

3. ローデータの可視化(Data Visualization)


フォロワー数の推移

フォロワー数の経時的な変化を可視化します。まずは、第一弾の特徴理解です。

# フォロワー数のトレンド
g1 <- ggplot(TidyAnime) + 
  theme_bw(base_size = 12, base_family = "HiraKakuProN-W6") +
  geom_line(aes(x = Time, y = Follower), colour = "darkgray") +
  geom_point(aes(x = Time, y = Follower)) +
  scale_x_datetime(expand = c(0.01, 0.01), breaks = scales::date_breaks("7 days")) +
  theme(axis.text = element_text(hjust = 1, angle = 45)) +
  scale_colour_manual(values = c("#6495ED", "#CD2626")) +
  labs(title = "ヴァイオレット・エヴァーガーデンのフォロワー数推移",
       subtitle = "放送開始1月10日-放送修了4月4日",
       x = "Date") +
  geom_vline(xintercept = c(as.POSIXct("2018-01-10"), as.POSIXct("2018-04-04")), linetype = 2)
# 描画
print(g1)