NTTドコモR&Dの技術ブログです。

【Python】YouTubeのトレンド分析:113ヵ国のデータを用いたトレンド傾向と生存時間分析(Cox比例ハザードモデル)

はじめに

本記事は、ドコモアドベントカレンダー2023の10日目の記事になります。

こんにちは。NTTドコモサービスイノベーション部4年目社員の九島です。

主な業務は、マーケティング分野への機械学習・AI活用です。

今回は、Kaggleで公開されているYouTubeのトレンド順位データ(113ヵ国分)を用いて、2部構成で分析を行いました。

  • 第1部では、YouTubeのトレンドデータから全体の傾向や国ごとの傾向、特に日本にフォーカスして分析してみました。
  • 第2部では、YouTubeのトレンドデータに生存時間分析を適用し、トレンドに上がった動画がトレンド外になるまでの持続時間を生存時間として分析してみました。

※本記事はPython環境が前提です。自分はJupyterノートブックで分析しました。

第1部:YouTubeトレンドデータの傾向分析

こちらの部では、YouTubeのトレンドデータから全体の傾向や国ごとの傾向、特に日本にフォーカスして分析してみました。

データセット

今回の分析には、Kaggleで公開されている以下のデータセットを使用しました。

Trending Youtube Video Statistics (113 Countries)

Kaggleサイトの右上のDownloadボタンからダウンロードし、自身の分析環境に配置します。

データの詳細は以下です。

カラム名 内容 データ型
title 動画のタイトル 文字列
channel_name チャンネル名 文字列
daily_rank トレンド順位 数値
daily_movement 前日からの変動順位 数値
weekly_movement 前週からの変動順位 数値
snapshot_date トレンド順位の取得日 文字列
country 国名コード 文字列
view_count 視聴回数 数値
like_count 高評価数 数値
comment_count コメント数 数値
description 説明文 文字列
thumbnail_url サムネイル画像のURL 文字列
video_id 動画のユニークID 文字列
channel_id チャンネルのユニークID 文字列
video_tags タグ 文字列
kind 種類 文字列
publish_date 動画投稿日 文字列
langauge 言語 文字列

試しに、日本のトレンド順位が1位の動画情報を可視化してみると、

import pandas as pd
data = pd.read_csv("./dataset/20231117/trending_yt_videos_113_countries.csv")
data[(data["country"]=="JP")&(data["daily_rank"]==1)].head()

私たちがよく存じている方々がお見えになりました!

このように、トレンド入りした動画の情報を確認することができます。

動画のトレンド順位は1~50位、それが113ヵ国それぞれに、取得された日数分あるというデータです。

また、今回使用するデータセットの一番古い日付が2023/10/26と割と最近でして(2023/12/10現在)、現時点ではあまり長くない期間のデータセットとなっています。

ただ、こちらのデータセットは毎日更新しているようでして(1日ずつデータが増えています)、1年後にまた分析してみたいなと思ったりもしますが、今回はこのまま分析します。

トレンド傾向の可視化

全体傾向

まずは全体の傾向から確認します。

データ内の数値データの統計量を見てみます。

data[["daily_movement","weekly_movement","view_count","like_count","comment_count"]].describe()

少しわかりずらいので、各カラムでヒストグラムを作成して見てみます。

import matplotlib.pyplot as plt
import seaborn as sns
%matplotlib inline
sns.set(font_scale=2)
plt.figure(figsize=(20, 10))
sns.distplot(data["daily_movement"], kde=False)
plt.show()

「daily_movement(前日からの変動順位)」は0付近に山があり、左右に割と綺麗に裾が広がっています。

マイナスの方が0に落ち込んでいることから、トレンド順位が急激に上がることは多く、急激に下がることは少ないと言えそうです(いきなりトレンド落ちはしていそうなので、50位以内での急激な低下は少ない、の方が表現として正しそうです)。

※ちなみに、初めてランクインした日のdaily_movementは、50-daily_rankで計算されていそうでして、その計算方法もこの分布に影響していそうです。

sns.set(font_scale=2)
plt.figure(figsize=(20, 10))
sns.distplot(data["weekly_movement"], kde=False)
plt.show()

「weekly_movement(前週からの変動順位)」はプラスの方向で一定に延びているような分布になっています。

これは一週間以上トレンド入りし続ける動画が少ないことを表していそうです。

※ちなみに、初めてランクインした日のweekly_movementも、50-daily_rankで計算されていそうでして、その計算方法もこの分布に影響していそうです。

sns.set(font_scale=2)
plt.figure(figsize=(20, 10))
sns.distplot(data["view_count"], kde=False)
plt.show()

こちらは「view_count(視聴回数)」の分布ですが、値のrangeが広く、このままだと見づらいので、対数を取ってみます。

import numpy as np
sns.set(font_scale=2)
plt.figure(figsize=(20, 10))
sns.distplot(np.log10(data["view_count"]+data[data["view_count"]!=0]["view_count"].min()), kde=False)
plt.show()

np.log10()をかけて、10の累乗で見てみますと、10の6乗付近で山になっていて、そこから左右に割と綺麗に裾が広がっています。

10の6乗は100万回の再生数なので、たしかにそれくらいの再生数は流行っている動画の感覚と一致しています。

sns.set(font_scale=2)
plt.figure(figsize=(20, 10))
sns.distplot(np.log10(data["like_count"]+data[data["like_count"]!=0]["like_count"].min()), kde=False)
plt.show()

こちらは「like_count(高評価数)」の分布で、先ほどのview_countと同様にlog10をかけています。

10の4乗から5乗の間で山になっています。

1万~10万なので、こちらも感覚とあっていそうです。

sns.set(font_scale=2)
plt.figure(figsize=(20, 10))
sns.distplot(np.log10(data["comment_count"]+data[data["comment_count"]!=0]["comment_count"].min()), kde=False)
plt.show()

最後に「comment_count(コメント数)」の分布で、先ほどと同様にlog10をかけています。

10の3乗付近で山になっています。

1000個もコメントあったら流行っている方ですね。

国ごとの傾向

続いて、国ごとの傾向を見ていきます。

groupby()を利用して同じ国のそれぞれの数値データの平均を取って、先ほどと同様にヒストグラムにしてみます(一気に図示します)。

sns.set(font_scale=2)
plt.figure(figsize=(20, 10))
sns.distplot(data.groupby("country").mean()["daily_movement"], kde=False ,bins=20)
plt.show()

sns.set(font_scale=2)
plt.figure(figsize=(20, 10))
sns.distplot(data.groupby("country").mean()["weekly_movement"], kde=False ,bins=20)
plt.show()

sns.set(font_scale=2)
plt.figure(figsize=(20, 10))
sns.distplot(data.groupby("country").mean()["view_count"], kde=False ,bins=20)
plt.show()

sns.set(font_scale=2)
plt.figure(figsize=(20, 10))
sns.distplot(data.groupby("country").mean()["like_count"], kde=False ,bins=20)
plt.show()

sns.set(font_scale=2)
plt.figure(figsize=(20, 10))
sns.distplot(data.groupby("country").mean()["comment_count"], kde=False ,bins=20)
plt.show()

どちらの数値もばらつきがあり、国ごとに異なる傾向になっているのは明らかでした。

また、それぞれ以下のような結果が確認できそうです。

  • daily_movement:2.5付近が多く、2~3の順位変化を平均とする国が多く、5や10を平均とする国が少しある。
  • weekly_movement:24の順位変化を平均とする国が非常に多く、18付近を平均とする国が少しある。
  • view_count:500万回以下を平均とする国が多く、1500万回付近を平均とする国が少しある。
  • like_count:20万以下を平均とする国が多く、20万~60万くらいを平均とする国が一定数ある。
  • comment_count:1万以下を平均する国が多く、1万~2万くらいを平均とする国が少しある。

さらに、各数値データのベスト1~3位とワースト1~3位の国を割り出してみました。

まず、ベスト1~3位を見てみます。

daily_movement = data.groupby("country").mean().sort_values("daily_movement",ascending=False)
weekly_movement = data.groupby("country").mean().sort_values("weekly_movement",ascending=False)
view_count = data.groupby("country").mean().sort_values("view_count",ascending=False)
like_count = data.groupby("country").mean().sort_values("like_count",ascending=False)
comment_count = data.groupby("country").mean().sort_values("comment_count",ascending=False)
quantitative_value_list = (daily_movement,weekly_movement,view_count,like_count,comment_count)
quantitative_value_name_list = ("daily_movement","weekly_movement","view_count","like_count","comment_count")
for q,q_name in zip(quantitative_value_list,quantitative_value_name_list):
    for number in range(3):
        print(f"{q_name}_best_{number+1}: ",q.index[number])

print出力結果は以下です。

daily_movement_best_1:  RU
daily_movement_best_2:  CA
daily_movement_best_3:  FR
weekly_movement_best_1:  CL
weekly_movement_best_2:  BR
weekly_movement_best_3:  GB
view_count_best_1:  KH
view_count_best_2:  MY
view_count_best_3:  LU
like_count_best_1:  MT
like_count_best_2:  LU
like_count_best_3:  MY
comment_count_best_1:  MT
comment_count_best_2:  EE
comment_count_best_3:  NZ

国名コードだとわかりにくいので、日本語に起こしてみます。

daily_movement_best_1:  ロシア連邦
daily_movement_best_2:  カナダ
daily_movement_best_3:  フランス共和国
weekly_movement_best_1:  チリ共和国
weekly_movement_best_2:  ブラジル連邦共和国
weekly_movement_best_3:  英国(グレートブリテン及び北アイルランド連合
王国)
view_count_best_1:  カンボジア王国
view_count_best_2:  マレーシア
view_count_best_3:  ルクセンブルク大公国
like_count_best_1:  マルタ共和国
like_count_best_2:  ルクセンブルク大公国
like_count_best_3:  マレーシア
comment_count_best_1:  マルタ共和国 
comment_count_best_2:  エストニア共和国 
comment_count_best_3:  ニュージーランド

意外なのか意外ではないのか判断が難しいですが、深堀りしがいのある結果になりました。

次に、ワースト1~3位を見てみます。

for q,q_name in zip(quantitative_value_list,quantitative_value_name_list):
    for number in range(3):
        print(f"{q_name}_worst_{number+1}: ",q.index[-(number+1)])

print出力結果は以下です。

daily_movement_worst_1:  HU
daily_movement_worst_2:  BG
daily_movement_worst_3:  LT
weekly_movement_worst_1:  LU
weekly_movement_worst_2:  ZW
weekly_movement_worst_3:  SV
view_count_worst_1:  BR
view_count_worst_2:  FR
view_count_worst_3:  DE
like_count_worst_1:  JP
like_count_worst_2:  BR
like_count_worst_3:  DE
comment_count_worst_1:  AL
comment_count_worst_2:  BA
comment_count_worst_3:  TZ

同様に日本語に起こしてみます(JPがありますね・・・)。

daily_movement_worst_1:  ハンガリー
daily_movement_worst_2:  ブルガリア共和国
daily_movement_worst_3:  リトアニア共和国
weekly_movement_worst_1:  ルクセンブルク大公国
weekly_movement_worst_2:  ジンバブエ共和国
weekly_movement_worst_3:  エルサルバドル共和国
view_count_worst_1:  ブラジル連邦共和国
view_count_worst_2:  フランス共和国
view_count_worst_3:  ドイツ連邦共和国
like_count_worst_1:  日本
like_count_worst_2:  ブラジル連邦共和国
like_count_worst_3:  ドイツ連邦共和国
comment_count_worst_1:  アルバニア共和国
comment_count_worst_2:  ボスニア・ヘルツェゴビナ
comment_count_worst_3:  タンザニア連合共和国

なんと、日本の高評価数の平均が最下位でした。

とても残念な結果ですが、ふと他の国を見てみると、ベストの~movementでランクインしていた国がワーストの~countでランクインしているのが確認できます。

この結果から、順位変動が大きい国、つまりトレンドの移り変わりが激しい国は、視聴回数や高評価数が少なくなる傾向にある可能性が浮上してきました。

試しに、日本の他の数値データの順位も見てみます。

print("daily_movement_JP: ",daily_movement.index.get_loc("JP")+1)
print("weekly_movement_JP: ",weekly_movement.index.get_loc("JP")+1)
print("view_count_JP: ",view_count.index.get_loc("JP")+1)
print("like_count_JP: ",like_count.index.get_loc("JP")+1)
print("comment_count_JP: ",comment_count.index.get_loc("JP")+1)

print出力結果は以下です(数値は順位を示しています)。

daily_movement_JP:  15
weekly_movement_JP:  8
view_count_JP:  110
like_count_JP:  113
comment_count_JP:  100

やはり、順位変動は大きく、回数系は少ない傾向がありました。

また、ベストとワーストの結果をよく見ると、この傾向は逆の場合も成り立ちそうです。

つまり、順位変動の大きさと回数系の多さは負の相関がありそうなことがわかりました。

たしかに、順位の変動が激しいと、次から次へと新しく面白い動画が投稿されていることになるので、視聴や高評価が分散してしまいそうですね。

日本でYouTube活動をするのは他の国と比べて大変なようです(あくまでトレンド入りさせるのを前提とした場合の話ですが)。

第2部:YouTubeトレンドの生存時間分析

こちらの部では、YouTubeのトレンドデータに生存時間分析を適用し、トレンドに上がった動画がトレンド外になるまでの持続時間を生存時間として分析してみました。

生存時間分析とは

生存時間分析は、時間が経過するにつれて特定のイベント(例:機械の故障、患者の死亡)が発生する確率を分析する統計的手法です。

この分析の目的は、様々な要因がイベント発生までの時間にどのように影響するかを理解し、予測することにあります。

生存時間分析は以下の特徴があります。

  • 打ち切りデータ(イベントが観測終了時点まで起こらないケース)を扱える。
  • 経過時間とイベント発生の関係を調べられる。
  • 様々な変数(年齢、性別、治療法など)がイベントのリスクにどのように影響するか分析できる。

応用例は以下です。

  • 医療分野での患者の生存期間予測
  • 機械の故障予測
  • 顧客の購買行動や離反予測

使用したモデル

今回は、Cox比例ハザードモデルを使用します。

Cox比例ハザードモデルは、生存時間分析で最も広く用いられるモデルの一つです。

このモデルは、リスク要因が時間と共にどのように影響するかをモデル化する際に使用されます。

Cox比例ハザードモデルは以下の特徴があります。

  • 比例ハザード仮定:各変数の効果が時間と共に変わらない(ハザード比が一定)と仮定。
  • 半パラメトリック:基準ハザード関数の形を仮定しないため、柔軟性が高い。
  • 多変量解析:複数の独立変数を同時に考慮。

本記事では詳細な数式は省略し、分析過程・結果にフォーカスします。

生存時間とイベント発生のデータ作成

生存時間分析に必要となる生存時間(トレンド持続日数)とイベント発生(トレンド落ち)をデータ化します。

まず、日付を扱いやすいように、データ型を変換します。

data_r = data.copy()
data_r["publish_date"] = pd.to_datetime(data_r["publish_date"]).dt.date
data_r["snapshot_date"] = pd.to_datetime(data_r["snapshot_date"]).dt.date

また、複数の国で同じ動画がランクインされることがあるので、video_idとcountryを結合してユニークなIDを作成します。

data_r["unique_id"] = data_r["video_id"] + "_" + data_r["country"]

以上の前処理を経て、トレンド持続日数を算出します。

video_survival_days = data_r.groupby("unique_id")["snapshot_date"].nunique().reset_index()
video_survival_days = video_survival_days.rename(columns={"snapshot_date": "survival_days"})
video_survival_days.head()

残りのトレンド落ちも算出します。 データ取得最新日時点でトレンドにまだ残っている動画は0、それ以外は1にします。

video_last_day = data_r.groupby("unique_id")["snapshot_date"].max().reset_index()
video_last_day = video_last_day.rename(columns={"snapshot_date": "last_day"})
video_last_day["not_ranked"] = 1
video_last_day.loc[video_last_day["last_day"]==data_r["snapshot_date"].max(),"not_ranked"] = 0
video_last_day.head()

ここで、簡単に生存時間の全体傾向と日本の傾向を見てみます。

sns.set(font_scale=2)
plt.figure(figsize=(20, 10))
sns.distplot(video_info["survival_days"], kde=False)
plt.show()

こちらは全体の傾向です。

大半の動画は5日以内にトレンド落ちしています。

sns.set(font_scale=2)
plt.figure(figsize=(20, 10))
sns.distplot(video_info.loc[video_info["country"]=="JP","survival_days"], kde=False)
plt.show()

こちらは日本の傾向です。

大半の動画は2,3日でトレンド落ちしています。

全体と比較して、生存日数は短めなようです。

Cox比例ハザードモデルによる生存時間予測モデルの構築

今回は、データセット内の数値データのみを利用して予測モデルを構築します。

まず、先ほど作成した生存時間データとイベント発生データに、数値データを結合します。

video_info = pd.merge(video_survival_days,video_last_day[["unique_id","not_ranked"]],on="unique_id",how="left")
video_first_day = data_r.groupby("unique_id")["snapshot_date"].min().reset_index()
video_first_day = pd.merge(video_first_day, data_r[["daily_rank","unique_id","snapshot_date","daily_movement", "weekly_movement", "view_count", "like_count", "comment_count"]],on=["unique_id","snapshot_date"],how="left")
video_info = pd.merge(video_info, video_first_day,on="unique_id",how="left")
video_info_r = video_info[[
    'survival_days', 
    'not_ranked', 
    'daily_rank', 
    'view_count', 
    'like_count', 
    'comment_count'
]]
video_info_r.head()

ここで注意ですが、今回使用するデータには、daily_movementとweekly_movementは含めませんでした。

なぜかというと、今回は生存時間を予測するタイムスパンを動画が初めてトレンド入りした日と置いたので、daily_movementとweekly_movementが全く意味のない情報になってしまうからです。

このタイムスパンを例えば、動画がトレンド入りしてから一週間後と置いた場合であれば、意味のある情報になると思いますが、今回はその分析はしていませんので、上記のようなデータ加工となっています。

続いて、Cox比例ハザードモデルを適用して予測モデルを構築します。

from lifelines import CoxPHFitter
cph = CoxPHFitter()
cph.fit(video_info_r, duration_col='survival_days', event_col='not_ranked')

こちらでモデル構築完了です。

では、モデルの詳細を見てみます。

cph.print_summary()

色々項目ありますが、注目すべきは、真ん中あたりのcoefカラムです。

こちらのcoefは生存時間への各変数の寄与度合いを表します。

daily_rankが0.03なので、トレンド入り時の順位は生存時間をプラスへ寄与すると読めます(別の言い方をすると、順位の値が大きいほど、トレンド落ちしやすい)。

また、他の変数が軒並み0.00なので、生存時間への寄与はほとんどないと読めます。

なんとも味気ない結果になりました。

つまり、トレンド入りしたときに高い順位であるほど、生存時間が長くなり、一方で視聴回数や高評価数はほとんど影響しないようです。

最後に、構築したモデルを用いて、動画の生存時間を予測してみます。

result = cph.predict_survival_function(video_info_r)
video0 = result.iloc[:, 2507] # トレンド入り時の順位が1位の動画
video1 = result.iloc[:, 9319] # トレンド入り時の順位が25位の動画
sns.set(font_scale=2)
plt.figure(figsize=(20, 10))
plt.plot(video0.index, video0)
plt.plot(video1.index, video1)

予測結果を比較するため、トレンド入り時の順位が1位(青色の線)と25位(橙色の線)の二つの動画の生存時間を予測してみました(日本の動画からピックアップしました)。

横軸が生存日数で、縦軸が生存確率を表しています。

結果から、トレンド入り時の順位が1位の動画が5日間生存する確率は50%で、25位の動画が5日間生存する確率は20%になり、トレンド入り時の順位が高いほど生存しやすいモデルになっていることが確認できました。

また、変数を変化させたときの生存時間の変化を可視化する機能もあるので、そちらも試してみます。

cph.plot_partial_effects_on_outcome(covariates='daily_rank', values=[5, 10, 15, 20, 25, 50], cmap='coolwarm')

これまでの結果と相違ないですね。

まとめ

本記事では、Kaggleで公開されているYouTubeのトレンド順位データ(113ヵ国分)を用いて、2部構成で分析を行いました。

  • 第1部では、YouTubeのトレンドデータから全体の傾向や国ごとの傾向、特に日本にフォーカスして分析してみました。
  • 第2部では、YouTubeのトレンドデータに生存時間分析を適用し、トレンドに上がった動画がトレンド外になるまでの持続時間を生存時間として分析してみました。

記事を読んでお気づきになられた方もいらっしゃると思いますが、今回の分析では文字列型のデータを使用していないため、少ない情報での分析に留まっておりました。

ですので、タイトル、説明文、サムネイル画像(URLから取ってくる)、動画ジャンル(説明文からなんとか抽出する)、チャンネル情報(チャンネルIDから引っ張ってくる)、等々のデータを分析に盛り込めれば、より深みのある分析が可能になると思っております。そちらは今後の課題とさせていただきたいと思います。

本記事が少しでもみなさまのお役に立てれば幸いです。

参考文献

Trending Youtube Video Statistics (113 Countries)

Pythonで生存時間解析 〜Cox比例ハザードモデルの解説と実装〜

Pythonで生存時間解析(Cox比例ハザードモデル, DeepSurv)