<p><div><div><div><div data-width="100" data-show-faces="false" data-send="true" data-layout="button_count" data-height="21" data-href="http://www.r-bloggers.com/recosystem-recommender-system-using-parallel-matrix-factorization/"></div></div><div><a data-url="http://www.r-bloggers.com/recosystem-recommender-system-using-parallel-matrix-factorization/" data-counturl="http://www.r-bloggers.com/recosystem-recommender-system-using-parallel-matrix-factorization/" data-count="horizontal" data-text="recosystem: Recommender System Using Parallel Matrix Factorization" class="twitter-share-button" data-via="rbloggers" href="https://twitter.com/share"></a></div><div></div></div><div></div></div>
<div>
(This article was first published on <strong><a href="http://statr.me/2016/07/recommender-system-using-parallel-matrix-factorization/"> Yixuan's Blog - R</a></strong>, and kindly contributed to <a rel="nofollow" href="http://www.r-bloggers.com/">R-bloggers)</a> </div>
<h1>A Quick View of Recommender System</h1>
<p>The main task of recommender system is to predict unknown entries in the<br>
rating matrix based on observed values, as is shown in the table below:</p>
<div align="center">
<img alt="Rating matrix" src="http://i0.wp.com/i.imgur.com/bmW79NS.png?w=456" data-recalc-dims="1">
</div>
<p>Each cell with number in it is the rating given by some user on a specific<br>
item, while those marked with question marks are unknown ratings that need<br>
to be predicted. In some other literatures, this problem may be named<br>
collaborative filtering, matrix completion, matrix recovery, etc.</p>
<p>A popular technique to solve the recommender system problem is the matrix<br>
factorization method. The idea is to approximate the whole rating matrix<br>
$R_{m\times n}$ by the product of two matrices of lower dimensions,<br>
$P_{k\times m}$ and $Q_{k\times n}$, such that</p>
<p></p>
<p>Let $p_u$ be the $u$-th column of $P$, and $q_v$ be the<br>
$v$-th column of $Q$, then the rating given by user $u$ on item $v$<br>
would be predicted as $p’_u q_v$.</p>
<p>A typical solution for $P$ and $Q$ is given by the following optimization<br>
problem [<a target="_blank" href="http://statr.me/2016/07/recommender-system-using-parallel-matrix-factorization/#FPSG2015">1</a>; <a target="_blank" href="http://statr.me/2016/07/recommender-system-using-parallel-matrix-factorization/#LRSG">2</a>]:</p>
<p></p>
<p>where $(u,v)$ are locations of observed entries in $R$, $r_{u,v}$ is<br>
the observed rating, $f$ is the loss function, and<br>
$\mu_P,\mu_Q,\lambda_P,\lambda_Q$ are penalty parameters<br>
to avoid overfitting.</p>
<p>The process of solving the matrices $P$ and $Q$ is called<br>
model training, and the selection of penalty parameters is<br>
parameter tuning. After obtaining $P$ and $Q$, we can then do<br>
the prediction of $\hat{R}_{u,v}=p’_u q_v$.</p>
<h1>LIBMF and recosystem</h1>
<p><a target="_blank" href="http://www.csie.ntu.edu.tw/~cjlin/libmf/">LIBMF</a><br>
is an open source C++ library for recommender system using parallel<br>
matrix factorization, developed by<br>
<a target="_blank" href="http://www.csie.ntu.edu.tw/~cjlin/">Dr. Chih-Jen Lin</a> and his research group.<br>
[<a target="_blank" href="http://statr.me/2016/07/recommender-system-using-parallel-matrix-factorization/#LIBMF">3</a>]</p>
<p>LIBMF is a parallelized library, meaning that<br>
users can take advantage of multi-core CPUs to speed up the computation.<br>
It also utilizes some advanced CPU features to further improve the performance.</p>
<p><a target="_blank" href="http://cran.r-project.org/package=recosystem">recosystem</a><br>
(<a target="_blank" href="https://github.com/yixuan/recosystem">Github</a>) is an R wrapper of<br>
the LIBMF library that inherits most of its features. Additionally, this<br>
package provides a number of user-friendly R functions to<br>
simplify data processing and model building. Also, unlike most other R packages<br>
for statistical modeling that store the whole dataset and model object in<br>
memory, LIBMF (and hence <code>recosystem</code>) can significantly reduce memory use,<br>
for instance the constructed model that contains information for prediction<br>
can be stored in the hard disk, and output result can also be directly<br>
written into a file rather than be kept in memory.</p>
<h1>Overview of recosystem</h1>
<p>The usage of <code>recosystem</code> is quite simple, mainly consisting of the following steps:</p>
<ol>
<li>Create a model object (a Reference Class object in R) by calling <code>Reco()</code>.</li>
<li>Specify the data source, either from a data file or from R objects in memory.</li>
<li>Train the model by calling the <code>$train()</code> method. A number of parameters<br>
can be set inside the function.</li>
<li>(Optionally) Call the <code>$tune()</code> method to select best tuning parameters<br>
along a set of candidate values, in order to achieve better model performance.</li>
<li>(Optionally) Export the model via <code>$output()</code>, i.e. write the factorization matrices<br>
$P$ and $Q$ into files or return them as R objects.</li>
<li>Use the <code>$predict()</code> method to compute predicted values.</li>
</ol>
<p>More details are covered in the package<br>
<a target="_blank" href="https://cran.r-project.org/web/packages/recosystem/vignettes/introduction.html">vignette</a><br>
and the help pages <code>?recosystem::Reco</code>, <code>?recosystem::data_source</code>, <code>?recosystem::train</code>,<br>
<code>?recosystem::tune</code>, <code>?recosystem::output</code>, and <code>?recosystem::predict</code>.</p>
<p>In the next section we will demonstrate how to use <code>recosystem</code> to analyze a<br>
real movie recommendation data set.</p>
<h1>MovieLens Data</h1>
<p>The <a target="_blank" href="https://movielens.org/">MovieLens</a> website collected many<br>
movie rating data for research use. [<a target="_blank" href="http://statr.me/2016/07/recommender-system-using-parallel-matrix-factorization/#MovieLens">4</a>] In this article we download the<br>
<a target="_blank" href="http://files.grouplens.org/datasets/movielens/ml-1m.zip">MovieLens 1M Dataset</a><br>
from <a target="_blank" href="http://grouplens.org/datasets/movielens/">grouplens</a>,<br>
which contains 1 million ratings from 6000 users and 4000 movies.</p>
<p>The rating data file, <code>ratings.dat</code>, looks like below:</p>
<div>
<pre><code>1::1193::5::978300760
1::661::3::978302109
1::914::3::978301968
1::3408::4::978300275
1::2355::5::978824291
...
</code></pre>
</div>
<p>Each line has the format <code>UserID::MovieID::Rating::Timestamp</code>,<br>
for example the first line says that User #1 gave Movie #1193 a rating of<br>
5 at certain time point.</p>
<p>In <code>recosystem</code>, we will not use the time information, and the required data<br>
format is <code>UserID MovieID Rating</code>, i.e., the columns are space-separated, and<br>
columns after <code>Rating</code> will be ignored.<br>
Therefore, we first transform the data file into a format that is supported by<br>
<code>recosystem</code>. On Unix-like OS’s, we can use the <code>sed</code> command to replace <code>::</code><br>
by a space:</p>
<div>
<pre><code>sed -e <span>'s/::/ /g'</span> ratings.dat > ratings2.dat
</code></pre>
</div>
<p>Then we can start to train a recommender, as the following code shows:</p>
<div>
<pre><code><span>library</span><span>(</span><span>recosystem</span><span>)</span><span> </span><span># 1
</span><span>r</span><span> </span><span>=</span><span> </span><span>Reco</span><span>()</span><span> </span><span># 2
</span><span>train_set</span><span> </span><span>=</span><span> </span><span>data_file</span><span>(</span><span>&quot;ratings2.dat&quot;</span><span>,</span><span> </span><span>index1</span><span> </span><span>=</span><span> </span><span>TRUE</span><span>)</span><span> </span><span># 3
</span><span>r</span><span>$</span><span>train</span><span>(</span><span>train_set</span><span>,</span><span> </span><span>opts</span><span> </span><span>=</span><span> </span><span>list</span><span>(</span><span>dim</span><span> </span><span>=</span><span> </span><span>20</span><span>,</span><span> </span><span># 4
</span><span> </span><span>costp_l1</span><span> </span><span>=</span><span> </span><span>0</span><span>,</span><span> </span><span>costp_l2</span><span> </span><span>=</span><span> </span><span>0.01</span><span>,</span><span> </span><span># 5
</span><span> </span><span>costq_l1</span><span> </span><span>=</span><span> </span><span>0</span><span>,</span><span> </span><span>costq_l2</span><span> </span><span>=</span><span> </span><span>0.01</span><span>,</span><span> </span><span># 6
</span><span> </span><span>niter</span><span> </span><span>=</span><span> </span><span>10</span><span>,</span><span> </span><span># 7
</span><span> </span><span>nthread</span><span> </span><span>=</span><span> </span><span>4</span><span>))</span><span> </span><span># 8
</span></code></pre>
</div>
<p>In the code above, line 2 creates a model object such that the training<br>
function <code>$train()</code> can be called from it. Line 3 specifies the data<br>
source – a data file on hard disk. Since in our data user ID and movie ID<br>
start from 1 rather than 0, we use the <code>index1 = TRUE</code> options in the function.</p>
<p>The data can also be read from memory, if the UserID, MovieID and Rating<br>
columns are stored as R vectors. Below shows an alternative way to provide<br>
the training set.</p>
<div>
<pre><code><span>dat</span><span> </span><span>=</span><span> </span><span>read.table</span><span>(</span><span>&quot;ratings2.dat&quot;</span><span>,</span><span> </span><span>sep</span><span> </span><span>=</span><span> </span><span>&quot; &quot;</span><span>,</span><span> </span><span>header</span><span> </span><span>=</span><span> </span><span>FALSE</span><span>,</span><span>
</span><span>colClasses</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>rep</span><span>(</span><span>&quot;integer&quot;</span><span>,</span><span> </span><span>3</span><span>),</span><span> </span><span>&quot;NULL&quot;</span><span>))</span><span>
</span><span>train_set</span><span> </span><span>=</span><span> </span><span>data_memory</span><span>(</span><span>user_index</span><span> </span><span>=</span><span> </span><span>dat</span><span>[,</span><span> </span><span>1</span><span>],</span><span>
</span><span>item_index</span><span> </span><span>=</span><span> </span><span>dat</span><span>[,</span><span> </span><span>2</span><span>],</span><span>
</span><span>rating</span><span> </span><span>=</span><span> </span><span>dat</span><span>[,</span><span> </span><span>3</span><span>],</span><span> </span><span>index1</span><span> </span><span>=</span><span> </span><span>TRUE</span><span>)</span><span>
</span></code></pre>
</div>
<p>Line 4 to line 6 set the relevant model parameters: $k, \mu_P,\mu_Q,\lambda_P$,<br>
and $\lambda_Q$, and Line 7 gives the number of iterations. Finally as I have<br>
mentioned previously, LIBMF is a parallelized library, so<br>
users can specify the number of threads that will be working<br>
simultaneously via the <code>nthread</code> parameter. However, when <code>nthread > 1</code>,<br>
the training result is <strong>NOT</strong> guaranteed to be reproducible, even if<br>
a random seed is set.</p>
<p>Now everything looks good, except one inadequacy: the setting of tuning<br>
parameters is ad-hoc, which may make the model sub-optimal.<br>
To tune these parameters, we can call the <code>$tune()</code> function to test<br>
a set of candidate values and use cross validation<br>
to evaluate their performance. Below shows this process:</p>
<div>
<pre><code><span>opts_tune</span><span> </span><span>=</span><span> </span><span>r</span><span>$</span><span>tune</span><span>(</span><span>train_set</span><span>,</span><span> </span><span># 9
</span><span> </span><span>opts</span><span> </span><span>=</span><span> </span><span>list</span><span>(</span><span>dim</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>10</span><span>,</span><span> </span><span>20</span><span>,</span><span> </span><span>30</span><span>),</span><span> </span><span># 10
</span><span> </span><span>costp_l2</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>0.01</span><span>,</span><span> </span><span>0.1</span><span>),</span><span> </span><span># 11
</span><span> </span><span>costq_l2</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>0.01</span><span>,</span><span> </span><span>0.1</span><span>),</span><span> </span><span># 12
</span><span> </span><span>costp_l1</span><span> </span><span>=</span><span> </span><span>0</span><span>,</span><span> </span><span># 13
</span><span> </span><span>costq_l1</span><span> </span><span>=</span><span> </span><span>0</span><span>,</span><span> </span><span># 14
</span><span> </span><span>lrate</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>0.01</span><span>,</span><span> </span><span>0.1</span><span>),</span><span> </span><span># 15
</span><span> </span><span>nthread</span><span> </span><span>=</span><span> </span><span>4</span><span>,</span><span> </span><span># 16
</span><span> </span><span>niter</span><span> </span><span>=</span><span> </span><span>10</span><span>,</span><span> </span><span># 17
</span><span> </span><span>verbose</span><span> </span><span>=</span><span> </span><span>TRUE</span><span>))</span><span> </span><span># 18
</span><span>r</span><span>$</span><span>train</span><span>(</span><span>train_set</span><span>,</span><span> </span><span>opts</span><span> </span><span>=</span><span> </span><span>c</span><span>(</span><span>opts_tune</span><span>$</span><span>min</span><span>,</span><span> </span><span># 19
</span><span> </span><span>niter</span><span> </span><span>=</span><span> </span><span>100</span><span>,</span><span> </span><span>nthread</span><span> </span><span>=</span><span> </span><span>4</span><span>))</span><span> </span><span># 20
</span></code></pre>
</div>
<p>The options in line 9 to line 15 are tuning parameters. The tuning function<br>
will evaluate each combination of them and calculate the associated<br>
cross-validated RMSE. The parameter set with the smallest RMSE will be contained<br>
in the returned value, which can then be passed to <code>$train()</code> (Line 19-20).</p>
<p>Finally, we can use the model object to do predictions. The code<br>
below shows how to predict ratings given by the first 20 users<br>
on the first 20 movies.</p>
<div>
<pre><code><span>user</span><span> </span><span>=</span><span> </span><span>1</span><span>:</span><span>20</span><span>
</span><span>movie</span><span> </span><span>=</span><span> </span><span>1</span><span>:</span><span>20</span><span>
</span><span>pred</span><span> </span><span>=</span><span> </span><span>expand.grid</span><span>(</span><span>user</span><span> </span><span>=</span><span> </span><span>user</span><span>,</span><span> </span><span>movie</span><span> </span><span>=</span><span> </span><span>movie</span><span>)</span><span>
</span><span>test_set</span><span> </span><span>=</span><span> </span><span>data_memory</span><span>(</span><span>pred</span><span>$</span><span>user</span><span>,</span><span> </span><span>pred</span><span>$</span><span>movie</span><span>,</span><span> </span><span>index1</span><span> </span><span>=</span><span> </span><span>TRUE</span><span>)</span><span>
</span><span>pred</span><span>$</span><span>rating</span><span> </span><span>=</span><span> </span><span>r</span><span>$</span><span>predict</span><span>(</span><span>test_set</span><span>,</span><span> </span><span>out_memory</span><span>())</span><span>
</span><span>library</span><span>(</span><span>ggplot2</span><span>)</span><span>
</span><span>ggplot</span><span>(</span><span>pred</span><span>,</span><span> </span><span>aes</span><span>(</span><span>x</span><span> </span><span>=</span><span> </span><span>movie</span><span>,</span><span> </span><span>y</span><span> </span><span>=</span><span> </span><span>user</span><span>,</span><span> </span><span>fill</span><span> </span><span>=</span><span> </span><span>rating</span><span>))</span><span> </span><span>+</span><span>
</span><span>geom_raster</span><span>()</span><span> </span><span>+</span><span>
</span><span>scale_fill_gradient</span><span>(</span><span>&quot;Rating&quot;</span><span>,</span><span> </span><span>low</span><span> </span><span>=</span><span> </span><span>&quot;#d6e685&quot;</span><span>,</span><span> </span><span>high</span><span> </span><span>=</span><span> </span><span>&quot;#1e6823&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>xlab</span><span>(</span><span>&quot;Movie ID&quot;</span><span>)</span><span> </span><span>+</span><span> </span><span>ylab</span><span>(</span><span>&quot;User ID&quot;</span><span>)</span><span> </span><span>+</span><span>
</span><span>coord_fixed</span><span>()</span><span> </span><span>+</span><span>
</span><span>theme_bw</span><span>(</span><span>base_size</span><span> </span><span>=</span><span> </span><span>22</span><span>)</span><span>
</span></code></pre>
</div>
<div align="center">
<img alt="Predicted ratings" src="http://i1.wp.com/i.imgur.com/nFGyyaO.png?w=456" data-recalc-dims="1">
</div>
<h1>Performance</h1>
<p>To make the best use of <code>recosystem</code>, the parallel computing option <code>nthread</code><br>
should be used in the training and tuning step. Also, LIBMF and recosystem can<br>
make use of some advanced CPU features to speed-up computation, if you<br>
compile the package from source and turn on some compiler options.</p>
<p>To build <code>recosystem</code>, one needs a C++ compiler that supports<br>
the C++11 standard. Then you can edit <code>src/Makevars</code> (<code>src/Makevars.win</code> for Windows system)<br>
according to the following guideline:</p>
<ul>
<li>The default <code>Makevars</code> provides generic options that should apply to most<br>
CPUs.</li>
<li>If your CPU supports SSE3<br>
(<a target="_blank" href="http://en.wikipedia.org/wiki/SSE3">a list of supported CPUs</a>), add</li>
</ul>
<div>
<pre><code>PKG_CPPFLAGS += -DUSESSE
PKG_CXXFLAGS += -msse3
</code></pre>
</div>
<ul>
<li>If not only SSE3 is supported but also AVX<br>
(<a target="_blank" href="http://en.wikipedia.org/wiki/Advanced_Vector_Extensions">a list of supported CPUs</a>), add</li>
</ul>
<div>
<pre><code>PKG_CPPFLAGS += -DUSEAVX
PKG_CXXFLAGS += -mavx
</code></pre>
</div>
<p>After editing the <code>Makevars</code> file, run <code>R CMD INSTALL recosystem</code> to install <code>recosystem</code>.</p>
<p>The plot below shows the effect of parallel computing and the compiler option on<br>
the performance of computation. The y axis is the elapsed time of the model tuning<br>
procedure in the previous example.</p>
<div align="center">
<img alt="Performance" src="http://i1.wp.com/i.imgur.com/GfgShWZ.png?w=456" data-recalc-dims="1">
</div>
<h2>References</h2>
<div></div>
<p>[1] Chin, Wei-Sheng, Yong Zhuang, Yu-Chin Juan, and Chih-Jen Lin. 2015a. <a target="_blank" href="http://www.csie.ntu.edu.tw/~cjlin/papers/libmf/libmf_journal.pdf"><em>A Fast Parallel Stochastic Gradient Method for Matrix Factorization in Shared Memory Systems</em></a>. ACM TIST.</p>
<div></div>
<p>[2] Chin, Wei-Sheng, Yong Zhuang, Yu-Chin Juan, and Chih-Jen Lin. 2015b. <a target="_blank" href="http://www.csie.ntu.edu.tw/~cjlin/papers/libmf/mf_adaptive_pakdd.pdf"><em>A Learning-Rate Schedule for Stochastic Gradient Methods to Matrix Factorization</em></a>. ACM TIST.</p>
<div></div>
<p>[3] Lin, Chih-Jen, Yu-Chin Juan, Yong Zhuang, and Wei-Sheng Chin. 2015. <a target="_blank" href="http://www.csie.ntu.edu.tw/~cjlin/libmf/"><em>LIBMF: A Matrix-Factorization Library for Recommender Systems</em></a>.</p>
<div></div>
<p>[4] F. Maxwell Harper and Joseph A. Konstan. 2015. <a target="_blank" href="http://dx.doi.org/10.1145/2827872"><em>The MovieLens Datasets: History<br>
and Context</em></a>. ACM Transactions on Interactive Intelligent Systems (TiiS) 5, 4,<br>
Article 19 (December 2015), 19 pages.</p>
<p><div><div><div><div data-width="55" data-show-faces="false" data-send="true" data-layout="box_count" data-height="62" data-href="http://www.r-bloggers.com/recosystem-recommender-system-using-parallel-matrix-factorization/"></div></div><div><a data-url="http://www.r-bloggers.com/recosystem-recommender-system-using-parallel-matrix-factorization/" data-counturl="http://www.r-bloggers.com/recosystem-recommender-system-using-parallel-matrix-factorization/" data-count="vertical" data-text="recosystem: Recommender System Using Parallel Matrix Factorization" class="twitter-share-button" data-via="rbloggers" href="https://twitter.com/share"></a></div><div></div></div><div></div></div>
<div>
<div>To <strong>leave a comment</strong> for the author, please follow the link and comment on their blog: <strong><a href="http://statr.me/2016/07/recommender-system-using-parallel-matrix-factorization/"> Yixuan's Blog - R</a></strong>.</div>
<hr>
<a rel="nofollow" href="http://www.r-bloggers.com/">R-bloggers.com</a> offers <strong><a rel="nofollow" href="http://feedburner.google.com/fb/a/mailverify?uri=RBloggers">daily e-mail updates</a></strong> about <a title="The R Project for Statistical Computing" rel="nofollow" href="http://www.r-project.org/">R</a> news and <a title="R tutorials" rel="nofollow" href="http://www.r-bloggers.com/?s=tutorial">tutorials</a> on topics such as: <a title="Data science" rel="nofollow" href="http://www.r-bloggers.com/?s=data%20science">Data science</a>, <a title="Big Data" rel="nofollow" href="http://www.r-bloggers.com/?s=Big%20Data">Big Data, <a title="R jobs" rel="nofollow" href="http://www.r-users.com/">R jobs</a>, visualization (<a title="ggplot and ggplot2 tutorials" rel="nofollow" href="http://www.r-bloggers.com/?s=ggplot2">ggplot2</a>, <a title="Boxplots using lattice and ggplot2 tutorials" rel="nofollow" href="http://www.r-bloggers.com/?s=boxplot">Boxplots</a>, <a title="Maps and gis" rel="nofollow" href="http://www.r-bloggers.com/?s=map">maps</a>, <a title="Animation in R" rel="nofollow" href="http://www.r-bloggers.com/?s=animation">animation</a>), programming (<a title="RStudio IDE for R" rel="nofollow" href="http://www.r-bloggers.com/?s=RStudio">RStudio</a>, <a title="Sweave and literate programming" rel="nofollow" href="http://www.r-bloggers.com/?s=sweave">Sweave</a>, <a title="LaTeX in R" rel="nofollow" href="http://www.r-bloggers.com/?s=LaTeX">LaTeX</a>, <a title="SQL and databases" rel="nofollow" href="http://www.r-bloggers.com/?s=SQL">SQL</a>, <a title="Eclipse IDE for R" rel="nofollow" href="http://www.r-bloggers.com/?s=eclipse">Eclipse</a>, <a title="git and github, Version Control System" rel="nofollow" href="http://www.r-bloggers.com/?s=git">git</a>, <a title="Large data in R using Hadoop" rel="nofollow" href="http://www.r-bloggers.com/?s=hadoop">hadoop</a>, <a title="Web Scraping of google, facebook, yahoo, twitter and more using R" rel="nofollow" href="http://www.r-bloggers.com/?s=Web+Scraping">Web Scraping</a>) statistics (<a title="Regressions and ANOVA analysis tutorials" rel="nofollow" href="http://www.r-bloggers.com/?s=regression">regression</a>, <a title="principal component analysis tutorial" rel="nofollow" href="http://www.r-bloggers.com/?s=PCA">PCA</a>, <a title="Time series" rel="nofollow" href="http://www.r-bloggers.com/?s=time+series">time series</a>, <a title="finance trading" rel="nofollow" href="http://www.r-bloggers.com/?s=trading">trading</a>) and more...
<img height="1" alt="" width="1" src="http://feeds.feedburner.com/~r/RBloggers/~4/1Yx-mIGos5A"></a></div>