u++の備忘録

AtCoder Beginner Contest 118をPythonで解く

AtCoder Beginner Contest 118」に出て、3完でした。

atcoder.jp

Cは実装に確証が持てておらず、かつコードも汚いので解説などで勉強します。Dは方針が立ちそうでしたが「ちょうどN本のマッチ棒を使って」という条件で分からなくなりました。

f:id:upura:20190216223457p:plain

A - B +/- A(100点)

  • 要件通りに実装
A, B = list(map(int, input().split()))

if (B % A == 0):
    print(A + B)
else:
    print(B - A)

B - Foods Loved by Everyone(200点)

  • 「好き」で挙げられたindexを一つのlistに格納
  • collectionsで集計
import collections

N, M = list(map(int, input().split()))

a = []
for i in range(N):
    tmp = list(map(int, input().split()))
    a += tmp[1:]

c = collections.Counter(a)
ans = 0
for val in c.values():
    if val == N:
        ans += 1

print(ans)

C - Monsters Battle Royale(300点)

  • 実装に確証がないまま通した
  • 最初に、一番小さい数字で、その他全てに可能な限り攻撃
  • その状態でソート
  • 小さい数字2種類を取り出し、互いに攻撃させ合った結果を答えとする
  • 数字が1種類しかない場合の処理も書く
N = int(input())
A = list(map(int, input().split()))

As = sorted(A)
As = [a%As[0] if a%As[0] != 0 else As[0] for a in As]
As = sorted(As)

def judge(tmp_one, tmp_two):
    flg = 1

    while flg:
        if tmp_two%tmp_one == 0:
            tmp_two = tmp_one
        else:
            tmp_two = tmp_two%tmp_one

        if tmp_one%tmp_two == 0:
            tmp_one = tmp_two
        else:
            tmp_one = tmp_one%tmp_two
        if tmp_one == tmp_two:
            flg = 0

    return tmp_one

tmp_one = As[0]
for i in range(1, N):
    if i == N - 1:
        tmp_two = As[i]
    elif As[i] != As[0]:
        tmp_two = As[i]
        break

ans = judge(tmp_one, tmp_two)
print(ans)

D - Match Matching(400点)

  • 桁数を増やすのを最優先で、マッチ棒の消費が小さい数字を積極的に作っていく?
  • 「ちょうど」使い切る組み合わせの探索が必要そうで、方針が分からなくなった
import numpy as np

N, M = list(map(int, input().split()))
A = list(map(int, input().split()))
A = sorted(A)

cost = [2, 5, 5, 4, 5, 6, 3, 7, 6]
candidate = []

for i in range(len(A)):
    candidate.append(cost[A[i]-1])

A = np.array(A)
candidate = np.array(candidate)

sort_index = candidate.argsort()
candidate = candidate[sort_index]
A = A[sort_index]

ans = ""
for _ in range(N//candidate[0]-1):
    ans += str(A[0])

rest = N%candidate[0] + candidate[0]
print(rest)
print(ans)

復習[追記: 20190217]

B

setにして積集合を取るのが簡単でした

note.nkmk.me


C

単純に最大公約数で良かったです。詳しくは公式解説を参照。

import fractions


N = int(input())
A = list(map(int, input().split()))

 
ans = fractions.gcd(A[0], A[1])

for i in range(2, N):
    ans = fractions.gcd(ans, A[i])

print(ans)

D

「ちょうどN本」の類いは、動的計画法(DP)の頻出問題のようですね。精進します。

qiita.com

初心者向けスポーツ分析チュートリアル「目標達成に導くデータ分析」 | Sports Analyst Meetup #1

2019年2月24日開催「Sports Analyst Meetup #1」の「初心者向けスポーツ分析チュートリアル」を担当します。本記事では、事前に発表資料を公開します。

connpass.com

「目標達成に導くデータ分析」と題して、スポーツアナリストの定義や分析の基礎的な話をまとめました。

大まかな流れは以下の通りです。

  1. スポーツアナリストの定義を紹介
  2. 定義で触れられていた「情報収集・分析力」「伝達力」について、事業会社のアナリストとしての解釈を付与
  3. 実際に身近なデータで分析していく実例を紹介

発表内ではプログラミング要素を割愛していますが、GitHubでコードを公開しています。

github.com

事前にご質問・ご要望があれば、twitterやブログのコメントなどでご連絡ください。

当日多くの皆さまとお会いできることを楽しみにしています。

「まんがで読破」をkindleの11円セールで30冊買って全部読んだ

f:id:upura:20190213114147p:plain

名作・名著を漫画化している「まんがで読破」シリーズが、kindleのセールで30冊が11円になっていました。30冊を購入し、1週間くらいかけて全部読んだので、備忘録として残しておきます。

ちなみに、この11円セールは明日2月14日までらしいです。

Twitterでの感想スレッド

以下のツイートから始まるスレッドに、30冊分の感想を1、2文程度で記載しています。これから読みたい人もいると思うので、意図的に抽象的な内容にしています。


お気に入り7選

あまりにも短い記事になりそうだったので、特に気に入った7冊を挙げておきます。

『般若心経』

説明が一般向けに簡略化された点も一部気になりましたが、全般にとても良くまとまっていました。

罪と罰

原作の文庫本も所持しているのですが、上下巻で分厚く尻込みしていました。「現代の預言書」と言われてる理由が垣間見えます。

君主論

主張だけ読むのではなく、この主張が生まれた背景も合わせて考えると一層面白いです。

ドグラ・マグラ

本当に頭おかしくなります。

『分析心理学 自我と無意識』

端的にまとめる感想は難しいのですが、本シリーズの中でも指折りでオススメです。

方法序説

「我思う故に我あり」のデカルトの著作です。名言が生まれた背景を簡易に説明しており、後半に振り返りのまとめも挿入され、分かりやすいです。

雇用・利子および貨幣の一般理論

本シリーズの一般的な分量の2倍に当たる400ページですが、分かりやすく一気に読み切れてしまいました。経済学を復習したくなりました。

AtCoder「みんなのプロコン2019」をPythonで解く

「みんなのプロコン2019」に出て、3完でした。

atcoder.jp

C問題のテストケースが一つだけWAになる理由が分からず、1時間以上を浪費しました。レートが下がって悲しいです。

f:id:upura:20190209223656p:plain

A - Anti-Adjacency(100点)

  • 偶数・奇数で場合分け
  • "Yes"ではないことに注意する
N, K = list(map(int, input().split()))

if N%2 == 0:
    a = N/2
else:
    a = int(N/2) + 1

if a >= K:
    print("YES")
else:
    print("NO")

B - Path(200点)

  • 大学で学んだグラフ理論を思い出した
  • 今回の設定では、2つのnodeから2つのedge、残りの2つのnodeから1つのedgeが出ていれば良い
  • 全てのinputを一つのlistに格納し、collectionsで条件を満たすか確認
import collections

a = []
for i in range(3):
    a += (list(map(int, input().split())))

c = collections.Counter(a)
ans = c.most_common()

flg = "NO"

if (ans[0][1]) == 2:
    if (ans[0][1]) == 2:
        if (ans[2][1]) == 1:
            if (ans[3][1]) == 1:
                flg = "YES"

print(flg)

C - When I hit my pocket...(400点)

  • 「1枚を2回連続で計2枚を増やす」「2回分を使って1円を経由し(B - A)枚を増やす」の比較
  • 前者が得な場合は、ずっと1枚を増やし続けるので答えは「K + 1」
  • 後者が得な場合は、可能な限り(B - A)枚を増やす
    • まずは1円に交換可能な時点まで進める
    • この時点で残りはrest = (K - A + 1)回
      • restが偶数なら、(B - A)枚増やせるのは(rest / 2)回
      • restが奇数なら、(B - A)枚増やせるのは((rest - 1) / 2)回、余った1回で1枚増やす
K, A, B = list(map(int, input().split()))

if (A + 2 >= B):
    print(K + 1)
elif (A >= K):
    print(K + 1)
else:
    rest = K - (A - 1)
    if (rest%2 == 0):
        print((B - A) * int(rest / 2) + A)
    else:
        print((B - A) * int(rest / 2) + A + 1)

一つだけWAになって時間を溶かした実装

  • pythonのfloatの桁落ちが原因だった
  • "(int(rest - 1) / 2)" のところが、正しくは "(int( (rest - 1) / 2))"
  • もっと楽なのは "/" の代わりに "//" 演算子


K, A, B = list(map(int, input().split()))
 
cnt = 1
 
if K < A + 1:
    cnt += K
    print(cnt)
    exit()
 
rest = K - (A - 1)
cnt += (A - 1)
 
if (B - A) > 2:
    if (rest%2 == 0):
        cnt += (B - A) * int(rest / 2)
    else:
        cnt += (B - A) * (int(rest - 1) / 2)
        cnt += 1
else:
    cnt += rest
 
print(int(cnt))

f:id:upura:20190209225933p:plain

pandasで協定世界時(UTC)を日本標準時(JST)に変換

はじめに

前回の記事で、iOS標準アプリ「ヘルスケア」からデータを書き出しcsvに変換する方法を紹介しました。

upura.hatenablog.com

ただし、このデータでは協定世界時(UTC)に準拠した日時情報が格納されています。日本標準時(JST)で過ごしている人が日次単位の集計を実施する際は、協定世界時日本標準時に変換するのが望ましいです。

本記事では、前回作成したcsvをpandasで読み込み、新たに日本標準時(JST)の列を作る方法をまとめました。

既存データ(協定世界時)を眺める

最初に、実際に日時情報が協定世界時で格納されていると確認してみます。

df_test = pd.read_csv('../data/input/df.csv')
df_test['datetime'] = pd.to_datetime(df_test['datetime'])

df_test['date'] = df_test['datetime'].dt.date
df_test['time'] = df_test['datetime'].dt.time

データは、普通の平日の1日分に絞り込んでおきます。また、その日のうちの累積の歩数を可視化することにします。

from datetime import date
df_test = df_test[df_test['date'] == date(2019, 1, 22)]
df_test['cumsumStepCount'] = df_test['stepCount'].cumsum()
plt.figure(figsize=(25, 12))
plt.rcParams["font.size"] = 18
plt.ylabel('cumsumStepCount')
plt.plot(df_test['time'], df_test['cumsumStepCount'])

f:id:upura:20190206114857p:plain

活動時刻が午前0時から午後2時あたりの範囲に収まっています。この日は、通常の平日と同様に出勤と退勤をしていたので、この時刻が日本標準時ではないと分かりました。

日本標準時に変換

まず、csvを読み込む段階で index_col に日時情報の列を指定しておきます。これは後述する日本標準時に変換するメソッドが、indexに対して実行しないとエラーが出る*1ためです。

df = pd.read_csv('../data/input/df.csv', index_col='datetime')
df.index = pd.to_datetime(df.index, utc=True)

f:id:upura:20190206115721p:plain

このindexを基に、日本標準時(JST)に変換した列を作成します。もちろんindexを上書きしても良いと思います。

df['datetimeja'] = df.index.tz_convert('Asia/Tokyo')

あとは、先ほどと同様の条件で可視化してみます。

df['date'] = df['datetimeja'].dt.date
df['time'] = df['datetimeja'].dt.time

df = df[df['date'] == date(2019, 1, 22)]
df['cumsumStepCount'] = df['stepCount'].cumsum()
plt.figure(figsize=(25, 12))
plt.rcParams["font.size"] = 18
plt.ylabel('cumsumStepCount')
plt.plot(df['time'], df['cumsumStepCount'])

f:id:upura:20190206120029p:plain

午前10時ごろに出勤、午後18時ごろに退社している様子が確認できました。この日は社外の勉強会に参加して、午後10時過ぎに帰路についたようです。

おわりに

本記事では、協定世界時に準拠した日時情報が格納されたcsvをpandasで読み込み、新たに日本標準時の列を作る方法をまとめました。

そもそも日次単位のデータをきちんと眺めないと、日時情報が協定世界時で格納されていると知らないままに分析してしまう事態にもなりかねません。改めて、データの仕様を把握する必要性を実感する事例でした。

iOS標準アプリ「ヘルスケア」からデータを書き出しcsvに変換

はじめに

本記事では、iOS標準アプリ「ヘルスケア」からデータを書き出し、csvに変換する方法をまとめます。

データの概要

ヘルスケアアプリはiOSに標準で搭載され、日常の歩数などが記録されています。自分に身近なデータなので、分析の仮説も立てやすく、データ分析の題材として便利かと思います。

データの取り出し方

手順は以下の通りです。

  1. ヘルスケアアプリからXMLファイルを書き出す
  2. XMLファイルをcsvファイルに変換する

ヘルスケアアプリからXMLファイルを書き出す

まずはヘルスケアアプリからデータを書き出します。この時点でcsv形式になっているPythonなどで扱いやすいのですが、XMLファイルでしか書き出すことはできません。

f:id:upura:20190203184634p:plain

まずは、カレンダーの画面から右上の「人物アイコン」をクリックします。

f:id:upura:20190203174420p:plain

この画面の下部にある「ヘルスケアデータを書き出す」をクリックすると、いくつかの形式でデータが書き出せるようになります。

f:id:upura:20190203174417p:plain

  • メッセージ
  • Gmail
  • LINE
  • Slack

など

f:id:upura:20190203174440p:plain

例えばSlackに送ると、次のようにzip形式で圧縮されたXMLファイルがメッセージとして投稿されます。

f:id:upura:20190203174434p:plain

zipファイルを展開すると、以下のようになっています。「書き出したデータ.xml」に、歩数などのデータが格納されています。

f:id:upura:20190203174448p:plain

中身をエディタで閲覧すると、以下のような形式になっています。

f:id:upura:20190203174423p:plain

XMLファイルをcsvファイルに変換する

ここでは、Pythonを用いてXMLファイルをcsvファイルに変換します。

import xml.etree.ElementTree as ET

tree = ET.parse('../data/Appleヘルスケア書き出し/書き出したデータ.xml')
root = tree.getroot()

上記の操作で、rootという変数にデータが格納されています。

f:id:upura:20190203191456p:plain

あとは要件に応じてデータを加工するだけです。

例えば歩数データが欲しい場合は、typeが 'HKQuantityTypeIdentifierStepCount' に等しいデータのみを取り出します。

import pandas as pd

datetime = []
counts = []

for child in root:
    data = child.attrib

    try:
        if data['type'] == 'HKQuantityTypeIdentifierStepCount':
            datetime.append(data['startDate'])
            counts.append(data['value'])
    except:
        pass

df = pd.DataFrame({
    'datetime': datetime,
    'stepCount': counts
})

df['stepCount'] = df['stepCount'].astype(int)

一定の日時ごとに、歩数のデータが記録されています。

f:id:upura:20190203191713p:plain

最後に、csv形式で保存しておきました。既にデータ自体は読み込めているので、csv形式でファイルを保存せず、そのまま分析を続けるのも一つの手だと思います。

df.to_csv('../data/input/df.csv', index=False)

分析例

ここでは、分析の一例として、書き出した生データから、私の日次の歩数を可視化してみます。

import pandas as pd

df = pd.read_csv('../data/input/df.csv', index_col='datetime')
df.index = pd.to_datetime(df.index, utc=True)
df['datetimeja'] = df.index.tz_convert('Asia/Tokyo')

データは、日付+時刻の形式で格納されています。ただし世界標準時で格納されているため、東京の時刻に変更しておきます。

upura.hatenablog.com

日次でまとめ上げるため、日付のみを取り出した列を作成しました。その後、その列をkeyにgroupbyすれば、日次の集計が可能です。

df['date'] = df['datetimeja'].dt.date
daySummary = df.groupby(['date']).sum()

f:id:upura:20190203192303p:plain

このDataFrameを折れ線グラフで可視化すると、次のようになりました。

import matplotlib.pyplot as plt
%matplotlib inline
plt.figure(figsize=(20, 12))
plt.rcParams["font.size"] = 18
plt.ylabel('stepCount')
plt.plot(daySummary)

f:id:upura:20190203192357p:plain

2017年10月ごろと2019年1月ごろに大きく跳ねているのは、それぞれ以下の取り組みが原因ですね。きちんと集計できているようです。

upura.hatenablog.com

upura.hatenablog.com

おわりに

本記事では、iOS標準アプリ「ヘルスケア」からデータを書き出し、csvに変換する方法をまとめました。冒頭でも書いた通り、自分に身近なデータなので、データ分析の題材として扱いやすいかと思っています。

「オレシカナイトデータコンペティションvol.2」で準優勝でした

はじめに

1月26日に、AbemaTVの実データを使ったデータ分析のコンテストに参加し、準優勝でした。本記事では、運営の確認を経て、問題ない範囲でコンテストの概要や私の解法を紹介します。

cyberagent.connpass.com

togetter.com

コンテストの概要

オレシカナイトとは

オレシカナイトとは、「AbemaTV」や「Ameba」をはじめとしたサイバーエージェントグループが運営するメディアのアドテクノロジー開発を行うエンジニアの横断組織「Cyberagent Publisher adTechnology Associaion」が主催する技術者向けの勉強会です。

今回はオレシカナイトの一企画として、AbemaTVの実データを用いて、広告視聴数の予測モデルの精度を競うコンテストが開催されました。同様の枠組みでの開催は二度目で、前回は学生のみを対象に実施されたとのことです。

今回の参加者は15人でした。大半が社会人で、学生は3人でした。

競技概要

過去2カ月分のデータを基に広告のインプレッション数を予測するモデルを構築し、将来1週間分の値を当てるコンテストでした。競技時間は10:30〜17:15。データは300MB程度のcsvファイルで提供され、クラウドサービスなどの利用も可能でした。

  • train: 543万レコード
  • test: 59万レコード

ビジネス的な意義

  • 出稿の段階で広告主とAbemaTV側でインプレッション数の値に関して契約を交わしている
    • 例:広告配信期間の2週間の間に、合計100万インプレッションを獲得する
  • この値を保証している関係上、予測の上で広告を配信し、合意した値に満たないと契約違反になる
  • 予測をあまりにも上回ると、実質的な値引きになってしまいビジネス上AbemaTV側が損をする
  • 上にも下にも外さない予測モデルが重要になる

既存の説明変数

既存の説明変数としては、以下のような情報が格納されていました。

  • チャンネルやCMの配信位置など広告枠に関する情報
  • CMの開始時刻などの広告配信タイミングに関する情報
  • インプレッション数(trainのみ)

広告枠に関する情報は匿名化されており、CMの開始時刻も順序を変えないままに一定期間ずらされているとのことでした。

評価指標

評価指標としては、AbemaTVが実際の運用の際に使っている指標を用いました。以下、便宜上「スコア」と表現します。

スコアは、以下の手順で算出されます。

  1. ある一定の塊でレコードをまとめ上げる
  2. 塊ごとに「実測値の合計/予測値の合計」の値を計算
  3. 値が0.9〜1.1に収まった塊の割合が「スコア」になる

例えば、以下の4レコードが一つの塊になっている場合を考えます。この時、(ユーザに非公開の)実測値がある中で、次の予測値を提出した状況です。

レコード 予測値 実測値
0000001 80 100
0000002 110 100
0000003 100 100
0000004 120 100

この時、「実測値の合計」は400で、「予測値の合計」は410になります。「実測値の合計/予測値の合計」は400/410=0.975となり、この塊はスコアに貢献します。

このような評価指標を用いているのは、実際のビジネスでの運用において「個々のレコード単位の予測精度」よりも「一定の塊単位での予測精度」の方が価値を持つからだそうです。詳細部分は都合上割愛しますが、広告主側・AbemaTV側が扱いやすい塊で集計が実施されているようです。

専用のオンライン採点ツールにcsvファイルを提出すると、スコアが表示され、リアルタイムで順位表が更新される仕組みとなっていました。予測値の提出は何度でも可能です。testデータの全ての値を用いたスコアが表示されていたので、不健全ではありますがtestに合わせたチューニングも可能な状況ではありました。

今回参加者側では、実際にどのような塊でtestのレコードがまとめ上げられているかは把握できませんでした。そこで参考値として、提出したファイルのスコアと合わせてMSEも表示されていました。

自分の解法

戦略

とにかく時間が短いコンテストだったので、Kaggle・Signateのコンペ、業務、ブログで過去に作成したコードを使い回す方針を徹底しました。

競技概要を聞いた際に、下記のSignateの公園コンペと類似していると感じました。実際に、特徴量作成のコードは、ほぼコピペで完了しています。適宜こちらの解説もご参照ください。

upura.hatenablog.com

端的に言えば「特徴量を作りLightGBMに投げ込む」という思考停止の戦略です。早い段階でベンチマークを作成し、中盤以降は探索的データ分析(EDA)もしながら、徐々に精度を高めていきました。

説明変数

  • 既存の広告枠に関する情報をカテゴリ変数として投入
  • CMの開始時刻からはmonth, day, day of the week, hour, minuteを抽出
  • day, day of the week, hour + minuteは、循環性を基に(sin, cos)にencoding
  • target mean encoding(5 fold)

学習モデル

  • LightGBM(Single model, 5 fold)
  • 参考値としてはMSEが表示されていましたが、大きな値に影響を受けづらいRMSEを最適化した方がスコアに寄与すると考え、RMSEをmetricに指定しました
  • MLPも試しましたが、カテゴリ変数が多いせいか、全く上手くいきませんでした

結果

惜しくも2位という結果になりました。開始1時間弱で最終的に4位となるスコアを出し、その後も1位を走っていただけに、終盤の16:30に抜かれてしまったのが非常に悔しいです。

1位の @hiding_koukyo さんの解法は下記をご参照ください。

koukyo1213.hatenablog.com

所感

1日で全日程を終える都合上、朝9:45集合という、私にとっては平日の定時よりも早い集合時間でした。

f:id:upura:20190129105924j:plain

運営の方々のサポートは、端的に言って素晴らしかったです。Slackでのやり取りは迅速でしたし、飲食物で言えば弁当も豪華でデザートの提供もありました。

2位の賞品としては、人形を頂きました。

私は今回の題材となったAbemaTVのプレミアム会員で、「相棒」をはじめとして頻繁にサービスを利用していました。残念ながら匿名化されている部分が多くドメイン知識を活かせる余地は少なかったですが、自分に馴染みの深いデータに触れられる貴重な体験でした。コンテストの評価指標もビジネス要件に紐付いており、独特な部分にどう立ち向かうか思考を凝らすことができ、非常に楽しかったです。

学びある、とても素晴らしい1日を送ることができました。参加者の皆さま、運営の皆さまに改めてお礼申し上げます。