ProbSpace論文コンペ振り返り（LB: 4位, Private: 7位）

この記事の概要

ProbSpaceにて開催された、論文引用数予測コンペの振り返り記事となります。
https://prob.space/competitions/citation_prediction
結果は、LB: 4位、Private: 7位と、Shake Downしてしまい、少し悔しい結果に終わってしまいました。
https://twitter.com/apolo_cr/status/1376203176366370818
CVとLBが結構綺麗に相関していたので、結構自信あったのですがね…。
とはいえ、終了2日前にLB: 3位に居れたことは大きな自信となりました！
ProbSpace初メダルが金メダルというのも嬉しかった！

そして何より、テーブルコンペながら、自然言語処理要素や時系列要素がある点や、目的変数の無いtrainデータが大量に与えられており、それをどのように活かすかが問われる点など、工夫しがいのある楽しいコンペでしたね！

それでは、コンペの概要及び、このコンペで私がどのような取り組みをしたかを、初心者の役に立つように、まとめていきたいと思います。

私の最終モデルのコードは、Githubに載せてありますので、もし良かったら参考にしてみてください！また、Twitterのフォロー及びZennのアカウントフォローもよろしくお願いします！（宣伝）
副業のお仕事依頼もお待ちしております。

https://github.com/apolo-cr/probspace_citation

コンペの概要

このコンペは、論文投稿サイト（プレ・プリントサーバー）の公開情報を用いて、被引用数を予測するというコンペでした。

上記はtrainデータで、citesが予想する目的変数なのですが、一番上がNaNとなっているのが分かりますね。このように、目的変数の無いtrainデータが大量あり、これをどう利用するかが課題となっていました。

また、与えられているデータとしては、title,abstract等のテキストデータやカテゴリ変数の他に、低精度被引用数doi_citesが与えられているのが特徴的なコンペでした。

Digital Object Identifier（DOI）により計算された低精度被引用数が代替変数として付与されています。(公式サイトより引用)

doi_citesをそのまま提出しても、LB: 0.7358のスコアが出るということで、随分と答えに近いデータが与えられていることが分かります。

評価指標は、RMLSEなので、（順当に？）目的変数のlogをとって、RMSEの最小化を目指すこととなります。

最終モデル概要

LightGBM×6、Catboost×3をRidge回帰でアンサンブルしたものが、Privateスコア最高（CV及びLBもほぼ最高に近い）でした！

CV Score	LB Score	Private Score
0.4827	0.4855	0.4902

作成した特徴量

特徴量は、groupbyのものや自然言語処理のもの全て含めて888個作成しました。

作成した特徴量のうち、特徴的なものを、いくつかピックアップして説明していきますね。

idとdoiの前半部分を擬似的なカテゴリ変数として利用

id列の「. 又は / の前」及び、doi列の「 / の前」を擬似的なカテゴリ変数として利用しました。どちらもコンペ取り組み初日から使っていたので、みんなやっていると思っていましたが、意外と使ってない人も多い(?)ようです。

追記：
5位のkatwooo414さんも同じことをしていたようですね！

idを「. または /」で区切ったPrefix（feature importanceはそれなりに高かったがリークしていただけかも）
5th Place Solution (katwooo414さんの解法)より引用

katwooo414さんはリークしていただけかもと言っていますが、リークはしていないと個人的には思っています。

ここで作成したカテゴリ変数は、以下のGroupby特徴量や自然言語処理の際に利用しています。
（以下、idの前半部分はid_feat、doiの前半部分はdoi_featと呼んでいます。）

Groupby特徴量

このGroupby特徴量では、doi_citesのカテゴリ変数ごとの最大値、平均、標準偏差を求めています。
カテゴリ変数には、NaNの多いjournal-refや report-noの他、上記のid_featやdoi_featも含めています。(当然Categories等も！)

def aggregation(input_df, group_keys, group_value, agg_methods):
    output_df = input_df[['id']]
    for agg_method in agg_methods:
        for col in group_keys:
            if callable(agg_method):
                agg_method_name = agg_method.__name__
            else:
                agg_method_name = agg_method
            if input_df.cites_known.min() == 1:
                new_col = f"{agg_method_name}_{group_value}_grpby_{col}"
            else:
                new_col = f"ALL_{agg_method_name}_{group_value}_grpby_{col}"
            new_df = (input_df[[col] + [group_value]].groupby(col)[[group_value]].transform(agg_method))
            new_df.columns = [new_col]
            output_df = pd.concat([output_df, new_df], axis=1)

    return output_df

GROUPBY_COLS = ['submitter', 'categories', 'license', 'id_feat', 'doi_feat', 'doi_feat_2']
def get_agg_doi_cites_features(input_df):
    group_keys = GROUPBY_COLS + ['authors', 'journal-ref', 'report-no']
    group_value = 'doi_cites'
    agg_methods = ['max', 'mean', 'std']
    output_df = aggregation(input_df, 
                            group_keys=group_keys,
                            group_value=group_value,
                            agg_methods=agg_methods)
    return output_df

時系列データを活かしたGroupby特徴量

以下の関数では、「そのカテゴリ内で、その論文が書かれる前に書かれた論文のdoi_citesの平均」を表しています。やっていることとしては、カテゴリ変数ごとに、作成日時でSortし、doi_citesを一つずつずらしたlag特徴量を作成し、lag特徴量のその時点までの平均を求めています。

def aggregation_CumFeat(input_df, group_keys, group_value, ascendings):
    output_df = input_df[['id']]
    for col in group_keys:
        for ascending in ascendings:
            if input_df.cites_known.min() == 1:
                new_col = f"CumFeat_{group_value}_grpby_{col}_{ascending}"
            else:
                new_col = f"ALL_CumFeat_{group_value}_grpby_{col}_{ascending}"
            input_df['lag'] = input_df.sort_values('first_created_unixtime', ascending=ascending).groupby([col])[[group_value]].shift(1)
            cum = input_df[[col] + ['lag', 'first_created_unixtime']].sort_values('first_created_unixtime').groupby(col).lag.agg(['cumsum', 'cumcount'])
            new_df = pd.DataFrame(cum['cumsum'] / cum['cumcount'])
            new_df.columns = [new_col]
            output_df = pd.concat([output_df, new_df], axis=1)

    return output_df

GROUPBY_COLS = ['submitter', 'categories', 'license', 'id_feat', 'doi_feat', 'doi_feat_2']
def get_agg_CumFeat_features(input_df):
    group_keys = GROUPBY_COLS + ['authors', 'journal-ref', 'report-no']
    group_value = 'doi_cites'
    ascendings = [True, False]
    output_df = aggregation_CumFeat(input_df, 
                            group_keys=group_keys,
                            group_value=group_value,
                            ascendings=ascendings)
    return output_df

ちなみに、この特徴量は、KaggleのRiiidコンペの公開notebookで皆が作成していた特徴量が元アイデアなのですが、時系列要素のあるテーブルコンペでは結構スコア改善に繋がる印象です。

本コンペ参加者が解法を考えるに当たって大いに参考にしていたであろうatmaCup#10でも、自分はこのような特徴量を作成して、スコア改善に寄与していました。

上記のコードと似たような発想で、(少し分かりにくいかもですが)単に平均を取るのではなく、時間との加重平均を取るイメージの特徴量も作成しました。

def aggregation_FirstTimeWeighted(input_df, group_keys, group_value, ascendings):
    output_df = input_df[['id']]
    for col in group_keys:
        for ascending in ascendings:
            if input_df.cites_known.min() == 1:
                new_col = f"FirstTimeWeighted_{group_value}_grpby_{col}_{ascending}"
            else:
                new_col = f"ALL_FirstTimeWeighted_{group_value}_grpby_{col}_{ascending}"
            input_df['lag'] = input_df.sort_values('first_created_unixtime', ascending=ascending).groupby(col)[[group_value]].shift(1)
            cum1 = input_df[[col] + ['lag', 'first_created_unixtime']].sort_values('first_created_unixtime').groupby(col).lag.agg(['cumsum'])
            cum2 = input_df[[col] + ['from_first_created', 'first_created_unixtime']].sort_values('first_created_unixtime').groupby(col).from_first_created.agg(['cumsum'])
            new_df = pd.DataFrame(cum1['cumsum'] / cum2['cumsum'])
            new_df.columns = [new_col]
            output_df = pd.concat([output_df, new_df], axis=1)

    return output_df

GROUPBY_COLS = ['submitter', 'categories', 'license', 'id_feat', 'doi_feat', 'doi_feat_2']
def get_agg_FirstTimeWeighted_features(input_df):
    group_keys = GROUPBY_COLS + ['authors', 'journal-ref', 'report-no']
    group_value = 'doi_cites'
    ascendings = [True, False]
    output_df = aggregation_FirstTimeWeighted(input_df, 
                            group_keys=group_keys,
                            group_value=group_value,
                            ascendings=ascendings)
    return output_df

こちらは、(賞金の出ない初心者向けコンペatmaCup#8を除けば)初のコンペ参加となった去年12月のKaggle Riiidコンペで、初コンペ初メダルを決定づける要因となった発想の特徴量です。
本コンペでもスコア改善に寄与しました。
このような特徴量は、GBTモデルが見つけるのが苦手な特徴だからこそスコア改善に寄与するのだと理解しています。

Word2Vecを利用してカテゴリ変数を特徴ベクトル化

atmaCup#10にて、アップロードされた、Araiさんの神ディスカッションをもとに作成しています。

本コンペとatmaCup#10は参加者が結構被っており、atmaCup#10終了後にLBが大きく動き出したのは、これを皆さんが取り入れたためだと個人的に思っています。

アイデア自体の説明は、上の神ディスカッションに任せるとして、本コンペでの私がやったことをまとめていきますね。

まず、Categories列を見ると、明らかに、複数のカテゴリがスペースで区切られて構成されているので、splitで分解しました。

df['categories_space_list'] = df['categories'].apply(lambda x: list(x.split()))

また、.前が同じで、.後が異なるカテゴリがあることに気付くので、.で分解したものも用意します。

df['categories_comma_list'] = df['categories'].apply(lambda x: list(re.split('[. ]', x)))

これらの、categories_space_listやcategories_comma_listを擬似的な文章とみなして、word2vecでベクトル表現に直した後、文章ごとに平均を取りました。

先ほど作成したid_featやdoi_featと結合したもの等にも同じ処理をしています。

model_size = {
    'id_doi': 16,
    'categories_space': 16,
    'categories_comma': 32,
    'id_doi_categories_space': 64,
}

n_iter = 100
w2v_dfs = []
for _df, _df_name in zip(
        [id_doi_df, categories_space_df, categories_comma_df, id_doi_categories_space_df],
        ['id_doi', 'categories_space', 'categories_comma', 'id_doi_categories_space']
    ):

    with timer(f"Creating w2v for {_df_name}"):
        # Word2Vecの学習
        w2v_model = word2vec.Word2Vec(_df['target'].values.tolist(),
                                    size=model_size[_df_name],
                                    min_count=1,
                                    window=100,
                                    workers=1,
                                    iter=n_iter)

    with timer(f"Getting document vector for {_df_name}"):
        # 各文章ごとにそれぞれの単語をベクトル表現に直し、平均をとって文章ベクトルにする
        sentence_vectors = _df['target'].progress_apply(lambda x: np.mean([w2v_model.wv[e] for e in x], axis=0))
        sentence_vectors = np.vstack([x for x in sentence_vectors])
        sentence_vector_df = pd.DataFrame(sentence_vectors,
                                        columns=[f"ALL_{_df_name}_w2v_{i}" for i in range(model_size[_df_name])])
        sentence_vector_df.index = _df['id']
        w2v_dfs.append(sentence_vector_df)

Categories, authors_parsedの処理

こちらは、上のWord2Vecを利用した特徴量と同様、複数のワードが同じセルに入っているCategories,authors_parsedから上手く特徴を作れないかと考えて作成したものです。

やっていることとしては、それぞれの単語に対して、CountEncodingや、上記のaggregation_CumFeat関数のような関数を適用し、複数の単語が含まれているセルは、各単語に与えられた数値を平均するというアイデアです。

上記のWord2Vecの特徴量を作成する前にこちらを作成していたので、スコア改善に大きく寄与しましたが、Word2Vec特徴量があっても尚、寄与するかどうかは定かでは無いです。

※もっと綺麗な書き方があるかもしれません。

def multiple_word_categories_features(input_df, target_cols):
    output_df = input_df[['id']]

    for target_col in target_cols:

        if input_df.cites_known.min() == 1:
            new_cols = [f'Multiple_{target_col}_count', f'Multiple_{target_col}_cumave', f'Multiple_{target_col}_anticipator']
        else:
            new_cols = [f'ALL_Multiple_{target_col}_count', f'ALL_Multiple_{target_col}_cumave', f'Multiple_{target_col}_anticipator']

        duplicates_containing_list = list(itertools.chain.from_iterable(list(input_df[target_col])))
        count_dic = collections.Counter(duplicates_containing_list)

        def categories_count_sum(l):
            sum = 0
            for i in l: 
                sum += count_dic[i]
            return sum

        keys = list(set(duplicates_containing_list))
        values = [0]*len(keys)
        doi_cites_cumsum_dic = dict(zip(keys, values))
        doi_cites_cumcount_dic = dict(zip(keys, values))

        keys = list(input_df.id)
        values = [0]*len(keys)
        id_cumave_dic = dict(zip(keys, values))
        id_anticipator_dic = dict(zip(keys, values))


        for target_list, doi_cites, id in tqdm(zip(input_df.sort_values('first_created_unixtime')[target_col], input_df.sort_values('first_created_unixtime')['doi_cites'], input_df.sort_values('first_created_unixtime')['id'])):
            cumave = 0
            anticipator = 0

            for i in target_list:
                doi_cites_cumsum_dic[i] += doi_cites
                doi_cites_cumcount_dic[i] += 1
                cumave += doi_cites_cumsum_dic[i]/doi_cites_cumcount_dic[i]
                anticipator += count_dic[i]/doi_cites_cumcount_dic[i]
            
            id_cumave_dic[id] = cumave
            id_anticipator_dic[id] = anticipator

        new_df = pd.DataFrame()
        new_df[new_cols[0]] = input_df[target_col].apply(categories_count_sum)
        new_df[new_cols[1]] = input_df['id'].map(id_cumave_dic)
        new_df[new_cols[2]] = input_df['id'].map(id_anticipator_dic)
        output_df = pd.concat([output_df, new_df], axis=1)

    return output_df

fastText及びSciBertの学習済みモデルをtitle,abstract,commentsに適用して、TruncatedSVDで次元圧縮

こちらは皆さんやっているかと思いますが、ホワイトリストに載っていたfastText及び申請で許可が出ていたSciBertの学習済みモデルを利用しています。

abstractは64次元、titleは32次元、commentsは16次元に圧縮しています。

model_en = load_model("cc.en.300.bin")

with timer("SVD abstract fastText"):
    X = df["abstract"].progress_apply(lambda x: model_en.get_sentence_vector(x.replace("\n", "")))
    X = np.stack(X.values)
    svd = TruncatedSVD(n_components=64, random_state=SEED)
    X = svd.fit_transform(X)
    X_df = pd.DataFrame(X, columns=[f"SVD_abstract_fastText_{i}" for i in range(X.shape[1])])
    X_df["id"] = df["id"]
    svd_fastText = svd_fastText.merge(X_df, on="id", how="left")

doi_citesをRidge回帰で置き換える

こちらは、どういう事かと言いますと、
doi_citesのみを説明変数として、Ridge回帰にかけ、その予測値を新たなdoi_citesとすることで、元々のdoi_citesがcitesと比べて大きく差がある部分を是正しようというものです。（結果的に外れ値処理に近いことになるのかな？）

上を見ればわかるように、実際のcitesは、doi_citesの値よりも大きい傾向にあります。
Ridge回帰をかけることで、doi_citesの値をcitesの値により近づけてからモデルを回しました。

※最終的にアンサンブルした9モデルのうち、4モデルでこの処理を行ってから学習させています。LB, Privateともにこの処理を行ったモデルの方が若干良いので、アプローチとしては悪くなかったでしょうか。

t = 15117
X_stage0 = np.log1p(df['doi_cites'])[:t].values.reshape(-1,1)
y_stage0 = np.log1p(df['cites'])[:t].values.reshape(-1,1)

alphas_list = [0.001, 0.005, 0.01, 0.05, 0.1, 0.5, 1, 5, 10, 100]
clf = RidgeCV(alphas=alphas_list, cv=5)
clf.fit(X_stage0, y_stage0)

oof = clf.predict(X_stage0)
pred = clf.predict(np.log1p(df['doi_cites'])[t:].values.reshape(-1,1))

# 参考
rmse = np.sqrt(mean_squared_error(y_stage0, oof))
print(rmse) # 0.6148

df['doi_cites'][:t] = oof.reshape(15117)
df['doi_cites'][t:] = pred.reshape(59084)

その他の特徴量

その他は、特筆すべき特徴量はありませんが、
・CountEncoding
・正規表現でpages及びfiguresの数値を抜き出し
・テキストカラムの文字数、単語数
等々作成しています。

ちなみに、シングルモデルで最もスコアの良かったCatboostのモデル(CV:0.4850,LB:0.4867,Private:0.4916)の、特徴量上位は以下のようになっています。

上位は、自然言語処理の特徴量で独占されていますね。ただ、自然言語処理の特徴量を一切作成する前のモデルでも、LightGBMシングルモデルで、LB: 0.4926, Private: 0.4967のスコアが出ているので、テーブルコンペの定番であるGroupby特徴量をしっかり作ることも大事かなと思います。

また、上のモデルのseedごとのスコアは以下の通りです。

Rondom Seed Averagingによってスコアが改善されているのが分かりますね！

ちなみに、Catboostは、6 × Seed Averagingで24時間近くモデル回すのにかかっているので、6Seedで限界でしたが、LightGBMでは、ちゃんと(?)、16 × Seed Averagingまでやっています。

上手くいかなかったこと

・Target Encoding
・feature importance上位の特徴量のみを用いて再学習（今回は作成特徴量888個と、そこまで無闇やたらに特徴量を増やしていないので、全てスコア改善に寄与していたのですかね）
・Xgboost

反省

上位陣の解法を見ると、つらなみdoi_citesとcitesの差分を目的変数にしているようでした。

図の右上部分見ると、citesが大きいデータを小さく予測してしまっていることがわかります。
citesが大きなデータは他のレンジのデータよりも少なく、予測が難しいのではと考えました。

そこで、弱ラベルであるdoi_citesがある程度citesに近いことを考えると、cites自体を予測するよりも、doi_citesがどのくらいcitesを正しく予測できているか(citesとdoi_citesの差分)を予測させた方が学習しやすいのではという仮説をたて、差分予測器を構築するに至りました。
6th Solution(ぐぐりらにっき)より引用

「citesが大きいデータを小さく予測している」ことには自分も気づいていて、差分を予測させるアイデアもあったのですが、どうせあまり変わらないかと決めつけて実験せず……。

思いつくのと、実際に実験するに至るまでには、やはり大きな差がありますね。

そして、LB: 4位からPrivate: 7位にShake Downしてしまったわけですが、何故でしょう…。
今回のデータ量が、trainデータ<<testデータであったので、CVとLB両方のスコアを注視する必要があると考えていましたが、幸い完璧に相関していたので、Shakeしないと思っていたのですが…。

Random Seed Averagingも大量にしていて、StratifiedKFoldも何十パターンもの平均となっていますし、最終モデルも、パラメータやfold数を色々変えた9モデルの平均ですし…。

下に全コードを載せておくので、「ここがShakeした原因だよ！」とか、
それ以外にも、「ここのコードの書き方汚い！」とかあったらアドバイス下さると嬉しいです。

もっと強くなりたい！！

Code

全過程のコードをGithubにupしましたので、ぜひ参考にしてみてください。
https://github.com/apolo-cr/probspace_citation