はじめに

こんにちは、事業会社で働いているデータサイエンティストです。

この記事では、次元数を自動判定し、かつ次元数が時間と共に変化することを許容するword2vecモデルを提案してみます。

早速モデルの説明に入ります！

モデル説明

モデルの数理的な部分はこちらを参照してください：

言葉で説明すると、このモデルは普通のword2vecとあまり変わらないです。word2vecのモデリングの際に魔法の数字として扱われがちな次元数をモデルのパラメータの一部として推定するところが特徴的です。かつ、その次元数は時間と共に変動します。

なぜ次元数を時間と共に変動するのを許容するのかというと、こちらの記事で示したように、時期によって急に現れたり消えたりする次元があり、この変動の可視化自体にビジネス的・アカデミア的な意義があるからです。

では早速モデルの詳細に入ります！

単語ベクトル

単語ベクトルの事前分布には、正規分布を設定することが多いですが、本当に正規分布で良いのか？という疑問に回答するため、正規分布から成すディリクレ過程で設定しています。

まずは全体のハイパーパラメーター$\alpha$をガンマ分布からサンプリングします。

$$
\alpha \sim Gamma(0.001, 0.001)
$$

次に、棒折り過程でディリクレ過程を構築します。無限にあるクラスターの中のクラスターpについてこのように必要な変数をサンプリングします：

$$
\pi_{p}\sim Beta(1, \alpha)
$$

$$
p_{p} = \pi_{p} \prod\limits_{l=1}^{p - 1} (1 - \pi_{l})
$$

また、クラスターpの中心ベクトルとそのばらつきは下記のようにサンプリングされます：

$$
P_{latent,p} \sim Normal(0,1)
$$

$$
P_{\sigma,p} \sim Gamma(0.001, 0.001)
$$

単語Sのベクトルをサンプリングする際は、まず上記の棒折り過程で構築した無限次元の確率ベクトルをパラメーターとするカテゴリ分布から、インデックスをサンプリングします：

$$
\eta_{S} \sim Categorical(p)
$$

次に、サンプリングされた$\eta_{s}$を基に、ベクトルをサンプリングします：

$$
ベクトル_{S} \sim Normal(P_{latent, \eta_{S}}, P_{\sigma, \eta_{S}})
$$

この一連のサンプリングのプロセス分布を$G$で表します。

次元数

次元数は自己回帰的な階層棒折り過程で定式化します。まず、最初の年の重要度ベクトルを棒折り過程でサンプリングします。具体的には、

$$
\beta \sim Gamma(0.001, 0.001)
$$

をサンプリングし、次に棒折り過程で$年重要度ベクトル_{1,d}$を求める

$$
\delta_{1,d}\sim Beta(1, \beta)
$$

$$
年重要度ベクトル_{1,d} = \delta_{1,d} \prod\limits_{l=1}^{d - 1} (1 - \delta_{1,l})
$$

これをd = 1からd = ∞まで繰り返します。

次に、経時的変化の度合いを示すパラメータをサンプリングします：

$$
\theta\sim Gamma(0.001, 0.001)
$$

そして、Y年の重要度ベクトルを、Y-1年の重要度ベクトルを「事前分布」とする階層棒おり過程からサンプリングします。

$$
\delta_{Y,d} \sim Beta\left( \theta \delta_{Y-1,d}, \theta \left(1 - \sum_{l=1}^{d} \delta_{Y-1,l} \right) \right)
$$

$$
年重要度ベクトル_{Y,d} = \delta_{Y,d} \prod\limits_{l=1}^{d - 1} (1 - \delta_{Y,l})
$$

これもd = 1からd = ∞まで繰り返します。

どうして階層化するのかというと、前の年の次元重要度のベクトルをこの年の年重要度ベクトルの事前分布にしないと、次元重要度のベクトル同士の重なりがなくなり、うまく経時的変化を捉えられません。

単語生成

単語の生成は古典的なword2vecのように、ネガティブサンプリングを採用します。要するに、本物のデータの偽物のデータを混入して、本物のデータを識別する方法を学習させます。

具体的には、単語iについて

$$
本物フラグ_{i}\sim Bernoulli(logit(\sum\limits_{a=1}^{\infty} (単語ベクトル_{i, a}*(\sum\limits_{c \in 周辺_{i}}^{}コンテキストベクトル_{単語_{c},a}) * 年重要度ベクトル_{年_{i}, a}))
$$

文脈とその単語が登場した年を共変量に、単語の出現を予測する無限次元ロジスティック回帰です。

モデル実装

モデル推定用のStanコードです：

dynamic_word2vec.stan

functions {
  vector stick_breaking(vector breaks){
    int length = size(breaks) + 1;
    vector[length] result;
    
    result[1] = breaks[1];
    real summed = result[1];
    for (d in 2:(length - 1)) {
      result[d] = (1 - summed) * breaks[d];
      summed += result[d];
    }
    result[length] = 1 - summed;
    
    return result;
  }
  
  real weighted_inner_product(
    vector vec_a, vector vec_b, vector weight
  ){
    vector[size(vec_a)] result_record;
    for (i in 1:size(vec_a)){
      result_record[i] = vec_a[i] * vec_b[i] * weight[i];
    }
    return sum(result_record);
  }
  
  real entropy(
    vector distribution, real threshold
  ){
    vector[size(distribution)] result;
    
    for (i in 1:size(distribution)){
      if (distribution[i] < threshold){
        result[i] = 0.0;
      } else {
        result[i] = distribution[i] * log(distribution[i]);
      }
    }
    
    return -sum(result);
  }
  
  real embedding_prior_lpmf(
    array[] int word_seq,
    int start, int end,
    
    int group_type,
    
    array[] vector word_embedding,
    
    vector group,
    array[] vector group_embedding, vector group_sigma
  ){
    vector[end - start + 1] lambda;
    int count = 1;
    for (w in start:end){
      vector[group_type] case_vector;
      for (g in 1:group_type){
        case_vector[g] = log(group[g]) + normal_lpdf(word_embedding[w] | group_embedding[g], group_sigma[g]);
      }
      lambda[count] = log_sum_exp(case_vector);
      count += 1;
    }
    return sum(lambda);
  }
  
  real partial_sum_lpmf(
    array[] int flag,
    int start, int end,
    
    array[] int year,
    array[] int word,
    array[] int word_lead_1,
    array[] int word_lag_1,
    array[] int word_lead_2,
    array[] int word_lag_2,
    array[] int word_lead_3,
    array[] int word_lag_3,
    
    array[] vector dimension_year,
    array[] vector word_embedding,
    array[] vector word_context
  ){
    vector[end - start + 1] lambda;
    int count = 1;
    for (i in start:end){
      lambda[count] = weighted_inner_product(
        word_embedding[word[i]],  
        word_context[word_lead_1[i]] + word_context[word_lag_1[i]] + 
        word_context[word_lead_2[i]] + word_context[word_lag_2[i]] + 
        word_context[word_lead_3[i]] + word_context[word_lag_3[i]], 
        dimension_year[year[i]]);
      count += 1;
    }
    return bernoulli_logit_lupmf(flag | lambda);
  }
}
data {
  int<lower=1> word_type;
  int<lower=1> year_type;
  int<lower=1> group_type;
  int<lower=1> dimension_type;
  
  array[word_type] int<lower=1,upper=word_type> word_seq;
  
  int<lower=1> N;
  array[N] int<lower=1,upper=year_type> year;
  array[N] int<lower=1,upper=word_type> word;
  array[N] int<lower=1,upper=word_type> word_lead_1;
  array[N] int<lower=1,upper=word_type> word_lag_1;
  array[N] int<lower=1,upper=word_type> word_lead_2;
  array[N] int<lower=1,upper=word_type> word_lag_2;
  array[N] int<lower=1,upper=word_type> word_lead_3;
  array[N] int<lower=1,upper=word_type> word_lag_3;
  array[N] int<lower=0,upper=1> flag;
  
  int<lower=0> val_N;
  array[val_N] int<lower=1,upper=year_type> val_year;
  array[val_N] int<lower=1,upper=word_type> val_word;
  array[val_N] int<lower=1,upper=word_type> val_word_lead_1;
  array[val_N] int<lower=1,upper=word_type> val_word_lag_1;
  array[val_N] int<lower=1,upper=word_type> val_word_lead_2;
  array[val_N] int<lower=1,upper=word_type> val_word_lag_2;
  array[val_N] int<lower=1,upper=word_type> val_word_lead_3;
  array[val_N] int<lower=1,upper=word_type> val_word_lag_3;
  array[val_N] int<lower=0,upper=1> val_flag;
}
parameters{
  real<lower=0> dimension_global_alpha;
  
  real<lower=0> dimension_across_year_alpha;
  array[year_type] vector<lower=0, upper=1>[dimension_type - 1] dimension_year_breaks;
  
  real<lower=0> group_alpha;                                       // ディリクレ過程の全体のパラメータ
  vector<lower=0, upper=1>[group_type - 1] group_breaks;  // ディリクレ過程のstick-breaking representationのためのパラメータ
  
  vector<lower=0>[group_type] group_sigma;
  array[group_type] vector[dimension_type] group_embedding;
  array[word_type] vector[dimension_type] word_embedding;
  array[word_type] vector[dimension_type] word_context;
}
transformed parameters {
  array[year_type] simplex[dimension_type] dimension_year;
  simplex[group_type] group;
  
  for (t in 1:year_type){
    dimension_year[t] = stick_breaking(dimension_year_breaks[t]);
  }
  
  group = stick_breaking(group_breaks);
}
model{
  dimension_global_alpha ~ gamma(0.001, 0.001);
  
  dimension_year_breaks[1] ~ beta(1, dimension_global_alpha);
  
  dimension_across_year_alpha ~ gamma(0.001, 0.001);
  
  for (t in 2:year_type){
    for (d in 1:(dimension_type - 1)){
      dimension_year_breaks[t, d] ~ beta(dimension_across_year_alpha * dimension_year[t - 1, d], dimension_across_year_alpha * (1 - sum(dimension_year[t - 1, 1:d])));
    }
  }
  
  group_alpha ~ gamma(0.001, 0.001);
  
  group_breaks ~ beta(1, group_alpha);
  
  group_sigma ~ gamma(0.001, 0.001);
  
  for (g in 1:group_type){
    group_embedding[g] ~ normal(0, 1);
  }
  
  int grainsize = 1;
  
  target += reduce_sum(
    embedding_prior_lupmf, word_seq, grainsize,
    
    group_type,
    
    word_embedding,
    
    group,
    group_embedding, group_sigma
  );
  
  target += reduce_sum(
    partial_sum_lupmf, flag, grainsize,
    
    year,
    word,
    word_lead_1,
    word_lag_1,
    word_lead_2,
    word_lag_2,
    word_lead_3,
    word_lag_3,
    
    dimension_year,
    word_embedding,
    word_context
  );
}
generated quantities {
  array[dimension_type] real drawn_G;
  array[word_type] vector[group_type] estimated_eta;
  vector[year_type] entropy_per_year;
  real F1_score;
  {
    int sampled_group = categorical_rng(group);
    drawn_G = normal_rng(group_embedding[sampled_group], group_sigma[sampled_group]);
    
    for (w in 1:word_type){
      vector[group_type] case_vector;
      for (g in 1:group_type){
        case_vector[g] = log(group[g]) + normal_lpdf(word_embedding[w] | group_embedding[g], group_sigma[g]);
      }
      estimated_eta[w] = softmax(case_vector);
    }
    
    for (i in 1:year_type){
      entropy_per_year[i] = entropy(dimension_year[i], 0.00001);
    }
    
    array[val_N] int predicted;
    int TP = 0;
    int FP = 0;
    int FN = 0;
    
    for (i in 1:val_N){
      predicted[i] = bernoulli_logit_rng(
        weighted_inner_product(
          word_embedding[val_word[i]],  
          word_context[val_word_lead_1[i]] + word_context[val_word_lag_1[i]] + 
          word_context[val_word_lead_2[i]] + word_context[val_word_lag_2[i]] + 
          word_context[val_word_lead_3[i]] + word_context[val_word_lag_3[i]], 
          dimension_year[val_year[i]])
      );
      if (val_flag[i] == 1 && predicted[i] == 1){
        TP += 1;
      }
      else if (val_flag[i] == 0 && predicted[i] == 1){
        FP += 1;
      }
      else if (val_flag[i] == 1 && predicted[i] == 0){
        FN += 1;
      }
    }
    F1_score = (2 * TP) * 1.0/(2 * TP + FP + FN);
  }
}

モデル推定

ここではモデル推定のための前処理とモデルコンパイルを実施します。

前処理の詳細はこちらを参照してください。

`%>%` <- magrittr::`%>%`
mecabbing <- function(text){
  text |>
    stringr::str_replace_all("[^一-龠ぁ-んーァ-ヶーa-zA-Z]", " ") |>
    RMeCab::RMeCabC(1) |>
    unlist() %>%
    .[which(names(.) == "名詞")] |>
    stringr::str_c(collapse = " ")
}

#　データのダウンロード
full_corp <- quanteda.corpora::download("data_corpus_foreignaffairscommittee") |>
  quanteda::corpus_subset(speaker != "")

# 肩書きを含む人名の抽出
capacity_full <- full_corp |>
  stringr::str_sub(1, 20) |>
  stringr::str_replace_all("\\s+.+|\n", "") |> # 冒頭の名前部分の取り出し
  stringr::str_replace( "^.+[参事|政府特別補佐人|内閣官房|会計検査院|最高裁判所長官代理者|主査|議員|副?大臣|副?議長|委員|参考人|分科員|公述人|]君(（.+）)?$", "\\1")

# 肩書きの抽出
capacity <- full_corp |>
  stringr::str_sub(1, 20) |>
  stringr::str_replace_all("\\s+.+|\n", "") |> # 冒頭の名前部分の取り出し
  stringr::str_replace( "^.+?(参事|政府特別補佐人|内閣官房|会計検査院|最高裁判所長官代理者|主査|議員|副?大臣|副?議長|委員|参考人|分科員|公述人|君(（.+）)?$)", "\\1") |>  # 冒頭の○から，名前部分までを消去
  stringr::str_replace("（.+）", "") |>
  stringr::str_replace("^○.+", "Other")

# 年の抽出
record_year <- quanteda::docvars(full_corp, "date") |>
  lubridate::year() |>
  as.numeric()

# 並列処理の準備
future::plan(future::multisession, workers = 16)

sub_corp_df <- tibble::tibble(
  capacity_full = capacity_full,
  capacity = capacity,
  record_year = record_year,
  # テキストデータ形式をquantedaからbase Rのcharacter型に変換する
  text = full_corp |>
    as.list() |>
    unlist()
) |>
  dplyr::filter(dplyr::between(record_year, 1981, 2017)) |>
  dplyr::filter(capacity %in% c("委員", "大臣", "副大臣")) |>
  # 発言の冒頭の名前と肩書きを削除
  dplyr::mutate(
    text = stringr::str_remove_all(text, stringr::fixed(capacity_full))
  ) |>
  dplyr::mutate(text_mecabbed = furrr::future_map_chr(text, ~ mecabbing(.x), .progress = TRUE)) |>
  dplyr::mutate(
    text_id = dplyr::row_number()
  )

sub_corp_use_tokens <- sub_corp_df |>
  dplyr::pull(text_mecabbed) |>
  quanteda::phrase() |>
  quanteda::tokens() |>
  quanteda::tokens_remove("の") |>
  quanteda::dfm() |>
  quanteda::dfm_trim(min_termfreq = 100, max_termfreq = 30000) |>
  colnames()

word_master <- tibble::tibble(
  word = sub_corp_use_tokens
) |>
  dplyr::arrange(word) |>
  dplyr::mutate(
    word_id = dplyr::row_number()
  )

year_master <- sub_corp_df |>
  dplyr::select(record_year) |>
  dplyr::distinct() |>
  dplyr::arrange(record_year) |>
  dplyr::mutate(
    year_id = dplyr::row_number()
  )

sub_corp_long_df_lag_correct <- sub_corp_df |>
  dplyr::select(text_id, record_year, text_mecabbed) |>
  dplyr::mutate(
    word = text_mecabbed |>
      quanteda::phrase() |>
      quanteda::tokens() |>
      quanteda::tokens_remove("の") |>
      as.list()
  ) |>
  tidyr::unnest_longer(word) |>
  dplyr::filter(word %in% word_master$word) |>
  dplyr::select(text_id, record_year, word) |>
  dplyr::left_join(word_master, by = "word") |>
  dplyr::left_join(year_master, by = "record_year") |>
  dplyr::group_by(text_id) |>
  dplyr::mutate(
    word_lag_1 = dplyr::lag(word_id, 1),
    word_lead_1 = dplyr::lead(word_id, 1),
    word_lag_2 = dplyr::lag(word_id, 2),
    word_lead_2 = dplyr::lead(word_id, 2),
    word_lag_3 = dplyr::lag(word_id, 3),
    word_lead_3 = dplyr::lead(word_id, 3)
  ) |>
  dplyr::ungroup() |>
  tidyr::drop_na()

set.seed(12345)

sub_corp_long_df_lag_correct_and_neg <- sub_corp_long_df_lag_correct |>
  dplyr::mutate(ans = 1) |>
  dplyr::bind_rows(
    sub_corp_long_df_lag_correct |>
      dplyr::mutate(word_id = sample(sub_corp_long_df_lag_correct$word_id, length(sub_corp_long_df_lag_correct$word_id), replace = FALSE)) |>
      dplyr::mutate(ans = 0)
  )

set.seed(12345)
val_id <- c(
  sample(which(sub_corp_long_df_lag_correct_and_neg$ans == 1), 5000, replace = FALSE),
  sample(which(sub_corp_long_df_lag_correct_and_neg$ans == 0), 5000, replace = FALSE)
)

data_list <- list(
  word_type = nrow(word_master),
  year_type = nrow(year_master),
  group_type = 10,
  dimension_type = 30,
  
  word_seq = 1:nrow(word_master),
  
  N = nrow(sub_corp_long_df_lag_correct_and_neg[-val_id,]),
  year = sub_corp_long_df_lag_correct_and_neg$year_id[-val_id],
  word = sub_corp_long_df_lag_correct_and_neg$word_id[-val_id],
  word_lead_1 = sub_corp_long_df_lag_correct_and_neg$word_lead_1[-val_id],
  word_lag_1 = sub_corp_long_df_lag_correct_and_neg$word_lag_1[-val_id],
  word_lead_2 = sub_corp_long_df_lag_correct_and_neg$word_lead_2[-val_id],
  word_lag_2 = sub_corp_long_df_lag_correct_and_neg$word_lag_2[-val_id],
  word_lead_3 = sub_corp_long_df_lag_correct_and_neg$word_lead_3[-val_id],
  word_lag_3 = sub_corp_long_df_lag_correct_and_neg$word_lag_3[-val_id],
  flag = sub_corp_long_df_lag_correct_and_neg$ans[-val_id],
  
  val_N = nrow(sub_corp_long_df_lag_correct_and_neg[val_id,]),
  val_year = sub_corp_long_df_lag_correct_and_neg$year_id[val_id],
  val_word = sub_corp_long_df_lag_correct_and_neg$word_id[val_id],
  val_word_lead_1 = sub_corp_long_df_lag_correct_and_neg$word_lead_1[val_id],
  val_word_lag_1 = sub_corp_long_df_lag_correct_and_neg$word_lag_1[val_id],
  val_word_lead_2 = sub_corp_long_df_lag_correct_and_neg$word_lead_2[val_id],
  val_word_lag_2 = sub_corp_long_df_lag_correct_and_neg$word_lag_2[val_id],
  val_word_lead_3 = sub_corp_long_df_lag_correct_and_neg$word_lead_3[val_id],
  val_word_lag_3 = sub_corp_long_df_lag_correct_and_neg$word_lag_3[val_id],
  val_flag = sub_corp_long_df_lag_correct_and_neg$ans[val_id]
)

m_dw2v_init <- cmdstanr::cmdstan_model("dynamic_word2vec.stan",
                                  cpp_options = list(
                                    stan_threads = TRUE
                                  )
)

では、モデル推定に行きます！

------------------------------------------------------------ 
EXPERIMENTAL ALGORITHM: 
  This procedure has not been thoroughly tested and may be unstable 
  or buggy. The interface is subject to change. 
------------------------------------------------------------ 
Gradient evaluation took 9.77816 seconds 
1000 transitions using 10 leapfrog steps per transition would take 97781.6 seconds. 
Adjust your expectations accordingly! 
Begin eta adaptation. 
Iteration:   1 / 250 [  0%]  (Adaptation) 
Iteration:  50 / 250 [ 20%]  (Adaptation) 
Iteration: 100 / 250 [ 40%]  (Adaptation) 
Iteration: 150 / 250 [ 60%]  (Adaptation) 
Iteration: 200 / 250 [ 80%]  (Adaptation) 
Iteration: 250 / 250 [100%]  (Adaptation) 
Success! Found best value [eta = 0.1]. 
Begin stochastic gradient ascent. 
  iter             ELBO   delta_ELBO_mean   delta_ELBO_med   notes  
   100     -4916474.274             1.000            1.000 
   200     -4475723.769             0.549            1.000 
   300     -4230342.443             0.385            0.098 
   400     -4022133.509             0.302            0.098 
   500     -3879118.970             0.249            0.058 
   600     -3785128.928             0.212            0.058 
   700     -3718109.838             0.184            0.052 
   800     -3668348.745             0.163            0.052 
   900     -3629199.934             0.146            0.037 
  1000     -3597975.928             0.132            0.037 
  1100     -3571687.102             0.121            0.025 
  1200     -3549102.281             0.111            0.025 
  1300     -3529624.091             0.103            0.018 
  1400     -3512299.355             0.096            0.018 
  1500     -3496721.876             0.090            0.014 
  1600     -3482390.685             0.085            0.014 
  1700     -3469253.548             0.080            0.011 
  1800     -3456956.962             0.076            0.011 
  1900     -3445386.102             0.072            0.009   MEDIAN ELBO CONVERGED 
Drawing a sample of size 1000 from the approximate posterior...  
COMPLETED. 
Finished in  28745.9 seconds.

最後に、結果を保存します：

m_dw2v_summary <- m_dw2v_estimate$summary()

推定結果の可視化

さて、予測精度から確認しましょう：

> m_dw2v_summary |>
     dplyr::filter(stringr::str_detect(variable, "F1"))
# A tibble: 1 × 7
  variable  mean median      sd     mad    q5   q95
  <chr>    <dbl>  <dbl>   <dbl>   <dbl> <dbl> <dbl>
1 F1_score 0.648  0.648 0.00431 0.00424 0.641 0.656

うーん、なんかいまいちかもしれないですね。期間が違うのでアレなんですが、前回の記事よりF1スコアが若干低いです。

まず、ディルクレ過程の分布を可視化したいですが、

> m_dw2v_summary |>
     dplyr::filter(stringr::str_detect(variable, "group\\["))
# A tibble: 10 × 7
   variable      mean   median       sd      mad        q5      q95
   <chr>        <dbl>    <dbl>    <dbl>    <dbl>     <dbl>    <dbl>
 1 group[1]  0.00106  0.00100  0.000401 0.000351 0.000531  0.00180 
 2 group[2]  0.993    0.994    0.00164  0.00146  0.990     0.996   
 3 group[3]  0.000550 0.000419 0.000441 0.000300 0.000120  0.00142 
 4 group[4]  0.000358 0.000285 0.000301 0.000219 0.0000690 0.000931
 5 group[5]  0.000356 0.000261 0.000317 0.000195 0.0000739 0.000976
 6 group[6]  0.000183 0.000104 0.000236 0.000102 0.0000167 0.000646
 7 group[7]  0.000220 0.000156 0.000212 0.000122 0.0000326 0.000627
 8 group[8]  0.000278 0.000204 0.000246 0.000156 0.0000527 0.000696
 9 group[9]  0.00303  0.00289  0.00110  0.000961 0.00150   0.00510 
10 group[10] 0.000476 0.000395 0.000344 0.000274 0.000111  0.00114

どうやら全ての単語ベクトルが同じ正規分布から抽出されているようです。なので、ディリクレ過程のような複雑なノンパラメトリックモデリングをしなくても、実は正規分布でよかったかもしれません。

さて、学習データ期間は1981年から2017年までで、日本の国内政治も国際政治も大きく変わりました、何回も何回も。だから、推定されたword2vecの次元もきっと色々変化しそうですね！

m_dw2v_summary |>
  dplyr::filter(stringr::str_detect(variable, "dimension_year\\[")) |>
  dplyr::mutate(
    id = variable |>
      purrr::map(
        \(x){
          stringr::str_split(x, "\\[|\\]|,")[[1]][2:3]
        }
      )
  ) |>
  tidyr::unnest_wider(id, names_sep = "_") |>
  dplyr::mutate(
    dplyr::across(
      dplyr::starts_with("id"),
      ~ as.integer(.)
    )
  ) |>
  dplyr::left_join(
    year_master, by = c("id_1" = "year_id")
  ) |>
  ggplot2::ggplot() + 
  ggplot2::geom_line(ggplot2::aes(x = record_year, y = mean, color = as.factor(id_2))) + 
  ggplot2::geom_ribbon(ggplot2::aes(x = record_year, ymin = q5, ymax = q95, fill = as.factor(id_2)), alpha = 0.3) +
  ggplot2::scale_x_continuous(breaks = seq(min(year_master$record_year), max(year_master$record_year), by = 2))

おーい！安定しすぎやろ！しかも11次元しか動いてないやないかい！ソ連が崩壊して、何回も政権交代したぞ！！！ふざけんな！推定時間を返せ！！と言いたくなると思いますが、次元の重要度が時間と共にあまり変化しないということ自体が興味深い発見です。いわゆる「国益」、ないしは「国益に関連する言説」が国際社会と国内社会の構造に関わらず安定していることを示しています。筆者の過去の単語ベクトル自体が時間と共に変動することを許容するモデルが過学習している可能性を示唆しています。

最後に、何個かの単語をピックアップして、コサイン類似度基準で最も類似している単語を確認しましょう。

まずは単語ベクトル行列を作ります：

word_embedding <- m_dw2v_summary |>
  dplyr::filter(stringr::str_detect(variable, "word_embedding\\[")) |>
  dplyr::pull(mean) |>
  matrix(ncol = 30)

次に、コサイン類似度の計算を実施します：

word_master |> 
  dplyr::filter(word == "フランス" | 
                  word == "韓国" | 
                  word == "中国" | 
                  word == "アメリカ" |
                  word == "イスラエル" |
                  word == "農業" | 
                  word == "正義" | 
                  word == "植民" | 
                  word == "テロ" |
                  word == "海賊" |
                  word == "円" |
                  word == "侵略" |
                  word == "支配" |
                  word == "安全" |
                  word == "冷戦" |
                  word == "国連" |
                  word == "ドル" |
                  word == "男性" |
                  word == "女性" |
                  word == "途上") |>
  dplyr::pull(word_id) |>
  purrr::map(
    \(id){
      word_embedding |>
        nrow() |>
        seq_len() |>
        purrr::map_dbl(
          \(i){
            (word_embedding[id,] %*% word_embedding[i,])/sqrt(sum(word_embedding[id,]^2) * sum(word_embedding[i,]^2))
          }
        ) |>
        tibble::tibble(
          word = word_master$word,
          cos_sim = _
        ) |>
        dplyr::arrange(-cos_sim) |>
        dplyr::select(word) |>
        dplyr::slice_head(n = 20)
    }
  ) |>
  dplyr::bind_cols() |>
  dplyr::glimpse()
New names:
• `word` -> `word...1`
• `word` -> `word...2`
• `word` -> `word...3`
• `word` -> `word...4`
• `word` -> `word...5`
• `word` -> `word...6`
• `word` -> `word...7`
• `word` -> `word...8`
• `word` -> `word...9`
• `word` -> `word...10`
• `word` -> `word...11`
• `word` -> `word...12`
• `word` -> `word...13`
• `word` -> `word...14`
• `word` -> `word...15`
• `word` -> `word...16`
• `word` -> `word...17`
• `word` -> `word...18`
• `word` -> `word...19`
• `word` -> `word...20`
Rows: 20
Columns: 20
$ word...1  <chr> "アメリカ", "ソ連", "側", "ソビエト", "接近", "友好国", "南", "台湾", "対ソ", "東側", "対", "思惑", "ロシア", "朝鮮民主主義人民共和国", "諸島…
$ word...2  <chr> "イスラエル", "イラン", "アラブ", "過激", "イスラム", "フジモリ", "パキスタン", "グループ", "ニカラグア", "ウクライナ", "ミッテラン", "干渉",…
$ word...3  <chr> "テロ", "予防", "セキュリティー", "闘い", "救援", "根絶", "テロリズム", "海賊", "イラク", "未然", "人命", "テロリスト", "後方", "戦闘", "特措…
$ word...4  <chr> "ドル", "円", "億", "倍", "倍増", "兆", "割", "平均", "弱", "トン", "分の", "合計", "当たり", "計上", "伸び", "比率", "率", "単位", "半分", "…
$ word...5  <chr> "フランス", "ドイツ", "ポーランド", "イギリス", "フィンランド", "アルゼンチン", "南アフリカ", "ブラジル", "マレーシア", "カナダ", "西独", "イ…
$ word...6  <chr> "中国", "韓国", "対", "接触", "台湾", "自動車", "朝鮮民主主義人民共和国", "摩擦", "側", "ソ連", "同士", "隣国", "ソビエト", "フィリピン", "コ…
$ word...7  <chr> "侵略", "戦争", "唯一", "被爆", "非難", "武装", "断固", "武力", "明白", "世論", "ポル", "絶対", "脅威", "疑惑", "原爆", "悪", "警告", "紛争",…
$ word...8  <chr> "円", "ドル", "億", "倍", "兆", "割", "倍増", "平均", "当たり", "分の", "比率", "弱", "半分", "計上", "額", "出資", "合計", "戸", "伸び", "単…
$ word...9  <chr> "冷戦", "対決", "第一歩", "軍拡", "融合", "崩壊", "核軍縮", "廃絶", "地帯", "終結", "湾岸", "北東", "核保有", "シナリオ", "ペレストロイカ", "…
$ word...10 <chr> "国連", "安保理", "憲章", "総会", "入り", "決議", "常任", "採択", "制定", "憲法", "総長", "創設", "テロリズム", "附帯", "発足", "敵国", "全会…
$ word...11 <chr> "女性", "職", "男子", "子供", "家庭", "結婚", "難民", "福祉", "教師", "名誉", "見舞い", "生徒", "婦人", "妻", "女子", "登録", "教育", "児童",…
$ word...12 <chr> "安全", "確保", "維持", "遂行", "平和", "保障", "目的", "利用", "向上", "抑止", "集団", "存立", "行使", "生命", "拘束", "排除", "行為", "恒久…
$ word...13 <chr> "支配", "勢力", "中立", "独占", "民族", "主義", "裏", "権力", "介入", "内政", "利益", "競争", "直接的", "自主", "助長", "潜在", "主権", "秩序…
$ word...14 <chr> "植民", "民衆", "中国人", "学生", "日系", "層", "高齢", "死傷", "若者", "孤児", "留学", "革命", "外国", "就労", "見舞い", "セルビア", "出身",…
$ word...15 <chr> "正義", "中立", "支持", "主義", "人権", "堂", "アパルトヘイト", "民族", "思想", "精神", "内政", "公正", "統一", "派", "パレスチナ", "ポル", "…
$ word...16 <chr> "海賊", "後方", "治安", "バグダッド", "活動", "部隊", "救援", "監視", "銃", "捜索", "テロ", "戦闘", "補給", "偵察", "出動", "小規模", "自衛隊…
$ word...17 <chr> "男性", "年齢", "就労", "永住", "結婚", "賃金", "介護", "高校", "師", "生徒", "看護", "負傷", "親", "学校", "高齢", "休暇", "配偶", "旅行", "…
$ word...18 <chr> "農業", "為替", "通貨", "ベース", "実績", "メリット", "市場", "文化", "比較", "格差", "需要", "界", "急激", "傾向", "産品", "工業", "先進", "…
$ word...19 <chr> "途上", "通貨", "先進", "産品", "相場", "債務", "援助", "発展", "食糧", "金融", "ウエート", "為替", "活性", "農業", "中南米", "新興", "木材",…
$ word...20 <chr> "韓国", "同士", "台湾", "接触", "先方", "朝鮮民主主義人民共和国", "ソ連", "折衝", "側", "中国", "ロシア", "東シナ海", "独", "牛肉", "豪州", "…

このように、ベクトルはきちんと単語の意味・共起関係を捉えられています。なので、このモデルの妥当性はある程度信用でき、日本の衆議院内の外交関連の言説はあまり質的に変化していないと理解しても問題なさそうです。

結論

いかがでしょうか？この分析自体からはあまりこれといった示唆は出せないですが、このモデルは問題なく時間と共に次元を選定できることが分かりましたので、今後はこのモデルでより言説に大きな変化がありそうなデータを分析します！

最後に、私たちと働いて、データサイエンスの力で社会を改善したい方はこちらをご確認ください：