廿TT

譬如水怙牛過窓櫺 頭角四蹄都過了 因甚麼尾巴過不得

ニューラルネット風ポアソン回帰でクリック数の多そうな検索クエリを見つける

今日の川柳

ニューラルネットポアソン回帰でしつこく遊んでいる。

検索クエリのスペースで区切られたひとかたまりの文字列をキーワードとよぶことにする。

検索クエリはキーワードの組み合わせで作られるので、クエリごとに生起したキーワードにフラグを立てていけば検索クエリを符号化できる。

例えば、「ポアソン+回帰」というクエリと「ワイブル+回帰」というクエリがあったら次の表のようにまとめられる。

ポアソン ワイブル 回帰
ポアソン+回帰 1 0 1
ワイブル+回帰 0 1 1

この 0 か 1 かの行列を説明変数にして、検索クエリのクリック数  Y_n を見積もるタスクを考える。

モデルはこう。

 Y_n \sim \mathrm{Poisson}(\exp(\mathrm{logit}^{-1}(X_n w_1)w_2'))

w_1w_2 の列数 K はニューラルネットにおける隠れ層のユニット数に対応する。

今回は K=4 とした。

確率的勾配降下法で最適化する。

f:id:abrahamcow:20190220003202p:plain

無事、収束しました。

クリック数の推定値と実測値をプロットする。

f:id:abrahamcow:20190220003236p:plain

当てはまりはいい感じ。

次に、学習したモデルで新しい検索クエリのクリック数を見積もってみる。

新しい検索クエリは行列Xにランダムに0か1を埋めることで生成する。

生成したクエリでクリック数の予測値が高いものをみれば、新たなコンテンツ作りのヒントが得られるかもしれない。

1000個クエリを生成して、現状クリックがないクエリのうち、予測クリック数が大きいトップ20を以下に図示する。

f:id:abrahamcow:20190220003830p:plain

ゴミですね。

以下にRのコードを貼ります。

library(searchConsoleR)
library(tidyverse)
library(Rcpp)
sourceCpp("~/Documents/NNPlogis.cpp")
scr_auth()
sc_websites <- list_websites()
scdata <- search_analytics(sc_websites[1,1],
                           startDate = "2019-02-01",endDate = "2019-02-07",
                           dimensions = c("query"))

query_list <-sapply(scdata$query,function(x)strsplit(x," "))
q_len <- sapply(query_list, length)
vocab <- unique(unlist(query_list))
n_vocab <- length(vocab)
q_co <- t(sapply(query_list, function(X)1*(vocab %in% X)))

Y <- matrix(scdata$clicks)
N <- nrow(Y)
D <- ncol(Y)
M <- ncol(q_co)
K <- 4
set.seed(4321)
w1 = matrix(rnorm(M*K), M, K)
w2 = matrix(rnorm(D*K), D, K)

out <- NNPlogis(Y,q_co,w1,w2,alpha = 0.01,max_iter = 2000)

plot(out$loglik,type = "l",xlab = "step",ylab = "log-likelihood")

Ybar <- exp(plogis(q_co%*%out$w1)%*%t(out$w2))

plot(Y,Ybar,xlab = "observed", ylab = "fitted")
abline(0,1,lty=2)

X <- matrix(0,1000,M)
set.seed(1)
X[cbind(rep(1:1000,each=2),sample.int(M,2*1000,replace = TRUE))] <- 1
pred <- drop(exp(plogis(X%*%out$w1)%*%t(out$w2)))
query <- apply(X==1,1,function(x)paste(vocab[which(x)],collapse = "+"))

dfout <- tibble(query=query,pred=pred) %>% 
  left_join(scdata) %>% 
  dplyr::filter(is.na(clicks)) %>% 
  dplyr::arrange(desc(pred)) %>% 
  top_n(pred,n=20)

ggplot(dfout,aes(x=reorder(query,pred),y=pred))+
  geom_col()+
  coord_flip()+
  theme_bw(20,"Osaka")+
  xlab("")+ylab("clicks (predicted)")

メインの関数は無駄にRcppで書いてます。

// [[Rcpp::depends(RcppArmadillo)]]
// [[Rcpp::depends(RcppProgress)]]
#include <RcppArmadillo.h>
#include <RcppArmadilloExtensions/sample.h>
#include <progress.hpp>
using namespace Rcpp;

arma::mat logis(const arma::mat x){
  return 1/(1+exp(-x));
}

// [[Rcpp::export]]
Rcpp::List NNPlogis(const arma::mat & Y,
                    const arma::mat & X,
                    const arma::mat & w1ini,
                    const arma::mat & w2ini,
                    const double & alpha,
                    const int & max_iter){
  int N = Y.n_rows;
  // int D = Y.n_cols;
  // int M = X.n_cols;
  // int K = w1ini.n_cols;
  arma::uvec pool = arma::linspace<arma::uvec>(0, N-1, N);
  arma::mat w1 = w1ini;
  arma::mat w2 = w2ini;
  arma::vec loglik(max_iter);
  Progress prog(max_iter);
  for(int i=0;i<max_iter;i++){
    arma::uvec I = Rcpp::RcppArmadillo::sample(pool,N,false);
    arma::mat sY = Y.rows(I);
    arma::mat sX = X.rows(I);
    double ll = 0;
    for(int n=0;n<N;n++){
      arma::mat tmpY = sY.row(n);
      arma::mat tmpX = sX.row(n);
      arma::mat z = logis(tmpX*w1);
      arma::mat fx = exp(z*w2.t());
      arma::mat delta1 = tmpY-fx;
      w2 += (alpha)*delta1.t()*z;
      w1 += (alpha)*(tmpX.t()*((delta1*w2)%(z%(1-z))));
      ll += sum(sum(tmpY % log(fx) -fx));
    }
    loglik[i] = ll;
    prog.increment();
  }
  return Rcpp::List::create(Rcpp::Named("w1")=w1,
                            _["w2"]=w2,
                            _["loglik"]=loglik);
}