【Rで自然言語処理】自然言語処理から見る米国映画の変遷 その壱

GWにたっぷり映画を見ました。映画ってその時代の空気を映し出しているよなぁ、と前々からなんとなく感じていたのですが、はっきりこれ、というデータがなくてむずがゆい感じがしていました。ならばいっそ、自分の手で分析してやろうじゃないか、というのが今回のモチベーションです。

※ものすごく長くなりそうなので二つに分割します。

データがなくちゃ始まらない

分析するといっても、データがなくちゃ始まらない。
必要なデータとしては、映画のあらすじ、公開年、制作国、ジャンルくらいあればいいでしょう。IMDbを漁ってみたのですが、ダウンロードできるデータがちょっと微妙な感じだったのでやめにして、こちらのデータを使うことにしました。

CMU Movie Summary Corpus

データベースのメンテナンスはどうやら2012年あたりで終わっているようですが、使う分には問題なさそうです。ちなみに2016年の映画のデータがほんのり入っていますが、もしかしたら公開前データを使っているのかもしれません。
上記サイトを見ていただくとわかるように、おおよそ必要なデータはそろっています。ファイルそのものはTSV形式なのでRで読み込むにも便利です。

f:id:wanko_sato:20180506181500p:plain

ダウンロードした圧縮ファイルを解凍するとこんな感じになっています。このうち、使うのはmovie.metadata.tsvとplot_summaries.txtの二つのファイルです。文字コードはいずれもUTF-8になっています。Rに取り込むのに面倒だったのでメモ帳で文字コードANSIに変更してしまいました。
文字コード変換がうまくいかなかったので。。。

というわけで、データが手に入ったので早速分析に取り掛かりましょう。

肝心なのは前処理だ

データセットが二つあります。ということは、最終的に分析に使うのに何らかの形でマージする、あるいはマージできるようにしておく必要があります。データ本体を見るとわかりますが、Wikipedia movie IDがマージキーにできそうです。それを念頭におきながら作業を進めていきます。

使用するライブラリ

###############################################################
# load library

library(dplyr)
library(textmineR)
library(text2vec)
library(tm)
library(topicmodels)
library(fastcluster)
library(tidytext)

今回使用するライブラリはこんな感じ。dplyrが活躍します。
fastclusterも便利です。このライブラリは距離行列を作らなくても、ベクトルの集合から直接階層的クラスタリングをしてくれる優れものなのです。なぜ優れものかというと、距離行列ってものすごくメモリを食うんですよね。たとえば10,000個のベクトルデータがあったとすると、距離行列は10,000×10,000の行列になってしまうわけで。メモリがすくないPCだと一発アウトなので、メモリが足りないときは重宝します。
※とはいえ、今はほとんど64bitマシンなので、そこまでメモリ制限を考える必要はないのですが。。。

ファイル読み込み

###############################################################
# data from 
# http://www.cs.cmu.edu/~ark/personas/

# set file-path
path <- file.path(getwd(),"movieAnalysis","MovieSummaries")

# read metadata
metadata <- read.delim(file=file.path(path,"movie.metadata_ansi.tsv"),
                                       header=F,
                                       sep="\t",
                                       stringsAsFactors=F)

# read summary data
summaryData <- read.delim(file=file.path(path,"plot_summaries_ansi.txt"),
                                              header=F,
                                              sep="\t",
                                              stringsAsFactors=F,
                                              quote="\n")

かなり時間がかかるのでのんびり待ってください。
ファイルパスは適当なところでどうぞ。
TSVファイルを読み込むのにはread.delimをつかってsep="\t"とすれば良いです。summaryの方はquote="\n"も必要です。データ間の区切りが改行(?)で扱われているようだったので。
で、データ読み込みができたので、下処理をしていきます。

整形してデータフレームにする

なんだかんだいってもやっぱりデータフレームが便利なわけで。

###############################################################
# create data

# create data-frame
dataDF <- lapply(seq(nrow(metadata)),function(x){
  ID <- metadata[x,1]
  title <- metadata[x,3]
  year <- ifelse(metadata[x,4]=="",NA,strsplit(metadata[x,4],split="-")[[1]][1])
  if(metadata[x,8]=="{}"){
    country <- NA
  }else{
    country <- gsub(pattern="\\{|}",
                    replacement="",
                    metadata[x,8]) %>%
      strsplit(split=",|: ") %>%
      unlist()
    country <- country[seq(from=2,to=length(country),by=2)]
    country <- paste(country,collapse=",")
  }
  if(metadata[x,9]=="{}"){
    genre <- NA
  }else{
    genre <- gsub(pattern="\\{|}",
                  replacement="",
                  metadata[x,9]) %>%
      strsplit(split=",|: ") %>%
      unlist()
    genre <- genre[seq(from=2,to=length(genre),by=2)]
    genre <- paste(genre,collapse=",")
  }
  summaryText <- ifelse(is.null(summaryData[summaryData$V1==ID,2]),NA,summaryData[summaryData$V1==ID,2])
  
  out <- data.frame(ID=ID,
                    TITLE=title,
                    YEAR=year,
                    COUNTRY=country,
                    GENRE=genre,
                    SUMMARY=summaryText)
  return(out)
})
dataDF <- do.call(rbind,dataDF) 

ちょっと長いですがやってることは単純です。movie.metadataはRDBの形になっているので必要な要素だけ抽出しています。summaryは前述したようにWikipedia movie IDでmetadataと紐づけられるので、合致するIDをもつmetadataに引っ付けていっています。
ここでポイントなのが国とジャンルの処理です。国とジャンルはちょっと面倒な入り方をしていて、

head(metadata[,8],30)
 [1] "{/m/09c7w0: United States of America}"  "{/m/09c7w0: United States of America}"  "{/m/05b4w: Norway}"                    
 [4] "{/m/07ssc: United Kingdom}"             "{/m/0345h: Germany}"                    "{/m/09c7w0: United States of America}" 
 [7] "{/m/0hzlz: South Africa}"               "{/m/09c7w0: United States of America}"  "{/m/0jgd: Argentina}"                  
[10] "{/m/07ssc: United Kingdom}"             "{/m/09c7w0: United States of America}"  "{/m/03_3d: Japan}"                     
[13] "{/m/09c7w0: United States of America}"  "{/m/07ssc: United Kingdom}"             "{/m/01znc_: Turkey}"                   
[16] "{/m/09c7w0: United States of America}"  "{}"                                     "{/m/09c7w0: United States of America}" 
[19] "{/m/03f2w: German Democratic Republic}" "{/m/09c7w0: United States of America}"  "{/m/09c7w0: United States of America}" 
[22] "{/m/09c7w0: United States of America}"  "{/m/09c7w0: United States of America}"  "{/m/0345h: Germany}"                   
[25] "{/m/05vz3zq: Soviet Union}"             "{/m/09c7w0: United States of America}"  "{/m/0f8l9c: France, /m/03rjj: Italy}"  
[28] "{/m/03rk0: India}"

まず全体が{ }でくくられており、国が複数の場合はカンマ区切り、かつ国の前にコロン区切りでコードが入っている、という構成になっているのです。ジャンルも同じ入り方をしています。このままだと取り扱いが非常に面倒です。とはいえ、データフレームにぶちこむわけですから、変数としては一つで持ちたい。そこでどうしたかというと、

  if(metadata[x,8]=="{}"){
    country <- NA
  }else{
    country <- gsub(pattern="\\{|}",
                    replacement="",
                    metadata[x,8]) %>%
      strsplit(split=",|: ") %>%
      unlist()
    country <- country[seq(from=2,to=length(country),by=2)]
    country <- paste(country,collapse=",")
  }

この部分になるわけです。

  1. gsubで全体をくくっている{ }をとりのぞき
  2. カンマまたはコロンでstrsplitし
  3. unlistしてベクトルに戻して
  4. 偶数番目のベクトル要素を抽出し
  5. カンマ区切りで結合

という処理をかましているのです。すると、

head(dataDF$COUNTRY,30)
 [1] United States of America   United States of America   Norway                     United Kingdom            
 [5] Germany                    United States of America   South Africa               United States of America  
 [9] Argentina                  United Kingdom             United States of America   Japan                     
[13] United States of America   United Kingdom             Turkey                     United States of America  
[17] <NA>                       United States of America   German Democratic Republic United States of America  
[21] United States of America   United States of America   United States of America   Germany                   
[25] Soviet Union               United States of America   France,Italy               India                     
[29] New Zealand                United States of America  

という形になってとても扱いやすくなります。

いらないものを消す

今回は米国の映画を対象としています。また、トピックモデルを使うのでsummaryがないものは対象としません。加えて、時系列分析もしてみたいので制作年が入っていないものも対象外とします。これらいらないものを消していくと、

dataDFreduce <- subset(dataDF,
                       !is.na(dataDF$YEAR) &
                       !is.na(dataDF$COUNTRY) &
                       !is.na(dataDF$GENRE) &
                       !is.na(dataDF$SUMMARY))
dataDFreduce <- dataDFreduce[grep(pattern="America",dataDFreduce$COUNTRY),]

このような感じになり、もともと55,000件ほどあったデータが13,000件ほどまで減りました。このデータフレームを中心にして分析を行っていきます。

分析その壱 トピックモデル

はい、お気に入りのトピックモデルLDAです。トピックモデルと合わせて階層的クラスタリングもやっています。この組み合わせがとても気に入っています。

wanko-sato.hatenablog.com

やっぱり前処理が必要だ

summaryをLDAでモデリングしてトピックの割り当て確率を算出、その類似度でもってクラスタリングする方法をとります。LDAでモデリングするにあたってはDTM(Document Term Matrix)を作る必要があるため、その前処理が必要になります。とはいえ、この前処理はほぼ定型なので、いつもの通り、な感じです。

###############################################################
# topic models

# create vector of summaries
allTexts <- as.character(dataDFreduce$SUMMARY)

# preprocess
sw <- c("i", "me", "my", "myself", "we", "our", "ours",
        "ourselves", "you", "your", "yours", tm:::stopwords("en"))
preText <- tolower(allTexts)
preText <- tm::removePunctuation(preText)
preText <- tm::removeWords(preText,sw)
preText <- tm::removeNumbers(preText)
preText <- tm::stemDocument(preText, language = "en")

まずsummaryのデータだけを抽出したベクトルを作ります。その際、どうやらdataDFreduce$SUMMARYがfactor型になっているっぽいのでcharacter型に変換します。stringasfactor=Fにしてるのに、なぜかこうなってしまうんですよね。気を付けていないと謎のエラーに悩まされます。
次に

  1. 全部を小文字に変換
  2. 句読点や記号など(punctuation)を削除
  3. StopWord(I,me等の頻出しすぎる単語)を削除
  4. 数字を削除
  5. ステミング(活用を削除して語幹のみにする)

の処理を加えています。このあたりは自然言語処理をするにあたっての定型処理です。

DTMをつくる

# tokenize(split into single words)
it <- itoken(preText,
             preprocess_function = tolower,
             tokenizer = word_tokenizer,
             ids = paste(dataDFreduce$YEAR,dataDFreduce$TITLE,sep=":"))

# delete stopwords and build vocablary
vocab <- create_vocabulary(it,
                           stopwords = sw)

# word vectorize
vectorizer <- vocab_vectorizer(vocab)

# create DTM
dtm <- create_dtm(it, vectorizer)

これも定型処理です。特に深く考える必要もないところです。

LDAのモデルを構築する

# modeling by package "topicmodels"
model <- LDA(dtm, control=list(seed=37464847,verbose=10),k = 20, method = "Gibbs")

"topicmodels"ライブラリのLDAを使っています。別の記事にも書きましたが、"topicmodels"のLDAモデルが個人的には扱いやすくて気に入っています。
ここではk=20としてトピック数を20にしています。あれこれ試した結果、だいたい20くらいで十分のようです。30までいくとごちゃごちゃしすぎますし、10だとうまく分離できない感じ。最終的にはこの20のトピックへの割り当て確率を使って類似度を計算してクラスタリングするので、特定のトピックにだけ注目する、ということは基本しません。なのでこれで良いのです。
また、seedを設定していますので、同じデータであれば同じ結果が出るはずです。

必要なデータを抽出する

termsDF <- get_terms(model,100)
topicProbability <- data.frame(model@gamma,dataDFreduce$YEAR,dataDFreduce$TITLE)
topicProbabilityFilter <- topicProbability[as.numeric(as.character(topicProbability$dataDFreduce.YEAR))>=1970,]
topicProbabilitySort <- topicProbabilityFilter[order(as.numeric(as.character(topicProbabilityFilter$dataDFreduce.YEAR))),]

LDAのモデルができたので各トピックの主要単語および各summaryのトピックへの割り当て確率を抽出します。また、合わせて1970年以降のデータのみ抽出するようにしています。
※今更だけれど最初から1970年以降のデータでトピックモデルを構築すりゃよかったんじゃないかと反省。

階層的クラスタリングで分割

さて、どうして各トピックへの割り当て確率でもって階層的クラスタリングをするかというと、

head(topicProbabilitySort)
            X1         X2         X3         X4         X5         X6         X7         X8         X9        X10        X11
55  0.05851064 0.03723404 0.23936170 0.02659574 0.02659574 0.02659574 0.02659574 0.02659574 0.02659574 0.05851064 0.02659574
58  0.09151786 0.05580357 0.45312500 0.01116071 0.01562500 0.02008929 0.03794643 0.01562500 0.03348214 0.01562500 0.03794643
399 0.08381503 0.03757225 0.07803468 0.04913295 0.02601156 0.05491329 0.02023121 0.03179191 0.05491329 0.03179191 0.04335260
434 0.03125000 0.09375000 0.05625000 0.04375000 0.03125000 0.09375000 0.03125000 0.08125000 0.03125000 0.03125000 0.03125000
460 0.03888889 0.06111111 0.03888889 0.02777778 0.02777778 0.11666667 0.06111111 0.02777778 0.06111111 0.02777778 0.03888889
482 0.03571429 0.17857143 0.02472527 0.01373626 0.04120879 0.07417582 0.01373626 0.11263736 0.08516484 0.01923077 0.03021978
           X12        X13        X14        X15        X16        X17        X18        X19        X20 dataDFreduce.YEAR
55  0.03723404 0.06914894 0.07978723 0.05851064 0.03723404 0.02659574 0.03723404 0.02659574 0.04787234              1970
58  0.01562500 0.02455357 0.03348214 0.02455357 0.01116071 0.01562500 0.01116071 0.01562500 0.06026786              1970
399 0.05491329 0.02601156 0.05491329 0.03179191 0.03179191 0.02023121 0.01445087 0.04335260 0.21098266              1970
434 0.05625000 0.04375000 0.08125000 0.04375000 0.05625000 0.03125000 0.05625000 0.04375000 0.03125000              1970
460 0.05000000 0.05000000 0.03888889 0.03888889 0.02777778 0.02777778 0.03888889 0.06111111 0.13888889              1970
482 0.02472527 0.01923077 0.01373626 0.01923077 0.10164835 0.10164835 0.05769231 0.01373626 0.01923077              1970
            dataDFreduce.TITLE
55              Carter’s Army
58                  Hell Boats
399                    Bigfoot
434            Promise at Dawn
460   ...tick...tick...tick...
482 Lovers and Other Strangers

このように、x1~x20までの変数に入っている値が各トピックへの割り当て確率になっており、そのベクトルが近ければ近いほど類似した文書である、ということが言えるからです。
例えば、一番目のトピック(ここではx1)の確率が高いからとその確率の高い文書だけを集めたとして、他のトピックへの割り当て確率が皆低い、とは限らないわけです。とすると、類似性としては微妙なわけで、文書の類似性を見るのであれば、ベクトル全体の類似性を見るべきではないか、と考えるのです。

さて、実行してみましょう。

###############################################################
# hierachical clustring

# clustering
hier <- fastcluster::hclust.vector(topicProbabilitySort[,1:20],method = "ward")
clustAssign <- cutree(hier,k=10)
topicProbabilityAssign <- cbind(topicProbabilitySort,clustAssign)

手法はWard法、距離行列にはユークリッド距離(fastclusterのデフォルト)を使います。これだとほら、距離行列を作るひと手間が省けるでしょ?
で、cutreeで適当なところで切って(今回はクラスタ数が10になるように)、それをデータフレームに引っ付けます。これでクラスタリングは完了。ここから可視化のためのデータを作っていきます。

可視化データ作成 箱ひげ図用

その前に、各クラスタにどれくらいの文書が割り当たったかを数えてみましょう。

table(topicProbabilityAssign$clustAssign)

   1    2    3    4    5    6    7    8    9   10 
 182 5990  670  391  452  335  418  251  307  176

はい、クラスタ2が異様に多いですね。おそらくクラスタ10個では分離しきれなかった微妙なやつがいっぱいあったのでしょう。それはクラスタ数を増やせば分割できそうな気もしますが、今回は無視しておきます。そこに手をかけるのが今回の目的ではないので。

# all clusters for boxplot by cluster
clusterDFlist <- split(topicProbabilityAssign,as.factor(topicProbabilityAssign$clustAssign))
plotClusterAll <- lapply(clusterDFlist,function(x){
  lapply(seq(20),function(y){
    x[,y]
  })
})

データフレームをクラスタ番号でsplitしてリストにします。そのリスト内のデータフレームから各トピックを一個ずつ抜き出してベクトルにし、リストにまとめています。で、

boxplot(plotClusterAll[[1]])

とすると、クラスタ1の各トピックへの割り当て確率が箱ひげ図で見れます。

f:id:wanko_sato:20180506192311p:plain

クラスタ1はトピック3への割り当てが特に高いようです。

f:id:wanko_sato:20180506192346p:plain

クラスタ2は明確な傾向がなさそうですね。

f:id:wanko_sato:20180506192406p:plain

クラスタ3はトピック15。

というように全部確認していくと、

クラスタNo トピックNo
1 3
2 なし
3 15
4 9
5 14
6 13
7 12
8 16
9 1
10 4

となっていることが分かりました。では、それぞれのクラスタがどんな映画なのか確認してみましょう。

クラスタの内容を確認する

クラスタ2は雑多に混ざっているので無視します。

The Odessa File, A Bridge Too Far, Commandoなど、戦争映画っぽいですね。

Beware! The Blob, Demon Seed, The Howlingなど、SF・ホラー系でしょうか。

The Sting, Young Sherlock Holmes, Pulp Fiction, Jackie Brownなど、ミステリー・クライム系っぽい感じ。

Luther, Malcolm X, American History Xなど、政治・社会系かな。

The Hidden, Alien Nation, Red Dragon・・・ん~これよくわからなかったんだけれど、なんだろう?

Halloween, Hell Night, The Texas Chainsaw Massacre Part 2など、スプラッタホラーですね。

Rocky II, Major League, Tysonなど、スポーツ物。

Jaws, Tentacles, Splash,The Bear・・・動物ものというかなんといいうか・・・

Excalibur, Willow, Dragonheartなど、ファンタジー物っぽいですね。

という感じで、はっきりこれとわかるものもあれば悩んでしまうようなものもありますが、おおむねきれいにクラスタリングできているように思います(主観)。

可視化データ作成 時系列用

さて、今回の最後に、各クラスタに分割された映画の数が時系列でどう変化していったかを確認します。

# count by year
countByYear <- lapply(clusterDFlist,function(x){
  out <- table(x$dataDFreduce.YEAR)
  out <- out[order(names(out))]
  out <- out[names(out)>=1970]
})

これで各クラスタの各年毎の映画数が集計できるので、それをmatplotで折れ線グラフにします。

映画数時系列推移

matplot(names(countByYear[[1]]),
        countByYearMat,
        col=seq(10),
        pch=1,
        type="o",
        lty=c(rep(1,8),rep(2,2)))
legend("topleft",legend=seq(10),
       col=seq(10),
       pch=1,
       lty=c(rep(1,8),rep(2,2)))

f:id:wanko_sato:20180506195642p:plain

赤い実線のクラスタ2が多すぎて他のクラスタがさっぱりわからないグラフになってしまいました。このままだとちょっとやな感じなので、

matplot(names(countByYear[[1]]),
        countByYearMat,
        ylim=c(0,40),
        col=seq(10),
        pch=1,
        type="o",
        lty=c(rep(1,8),rep(2,2)))
legend("topleft",legend=seq(10),
       col=seq(10),
       pch=1,
       lty=c(rep(1,8),rep(2,2)))

※ylimを加えただけです。

f:id:wanko_sato:20180506195834p:plain

ちょっとごしゃごしゃしていますが、2002~4年あたりから急激に増えているクラスタが見えます。わかりやすいように取り出してみると、

matplot(names(countByYear[[1]]),
        countByYearMat[,c(1,3,5,7)],
        ylim=c(0,40),
        pch=1,
        type="o",
        col=c(1,3,5,7),
        lty=1)
legend("topleft",legend=c(1,3,5,7),
       col=c(1,3,5,7),
       pch=1,
       lty=1)

f:id:wanko_sato:20180506200014p:plain

クラスタ3,5,7が大きく伸びていて、一方でクラスタ1がさっぱり増えていません。クラスタ1の戦争ものは1980年代半ばくらいまではほかの系統と同じくらい作られていましたが、その後大きく数は増えていません。一方、クラスタ3のSF・ホラーは80年代半ばから大きく数を増やしており、2000年に入って間もなく、さらに数を増やしています。同じくホラーでもクラスタ7のスプラッタは90年代に入って数を減らしていましたが、2000年に入って盛り返してきました。クラスタ5の政治・社会ものが2000年に入って急激に増加したのは世紀が変わって社会的・政治的な変動が大きかった影響のように感じます。

こうしてみてくると、やっぱり映画って時代の雰囲気やらそのときの社会情勢やらを反映しているのだなぁ、としみじみ感じるのです。

まとめ

とりあえず今回はトピックモデルとクラスタリングでどんな内容の映画が作られているのか、それが時系列でどのように変化しているか、を見てきました。これだけだとちょっと踏み込みが浅い感じがしたので、さらにSentiment Analysisもやっています。が、今回はさすがに長くなりすぎたので分割して次回に回します。