From 23d501303a692bcd1b75b23739e6117f617e1d1a Mon Sep 17 00:00:00 2001 From: Norm Matloff Date: Wed, 10 May 2017 21:00:03 -0700 Subject: [PATCH 01/30] improved doc --- README.md | 34 ++++++++++++++++++---------------- 1 file changed, 18 insertions(+), 16 deletions(-) diff --git a/README.md b/README.md index ba45ef3..a106b29 100644 --- a/README.md +++ b/README.md @@ -9,7 +9,8 @@ * Parallel computation. -* Novel variations on statistical latent factor model. +* Novel variations on common models, e.g. a hybrid of + NMF and k-Nearest Neighbor. * Plotting. @@ -27,9 +28,9 @@ A simple random effects latent factor model is E(Y) = μ + αi + βj -where Yij is the rating, α and βj being -specific latent effects for user i and item j, e.g. movie reviewer i and -movie j. +where Yij is the rating, with αi and +βj being specific latent effects for user i and item j, +e.g. movie reviewer i and movie j. Though typically Maximum Likelihood Estimation is used for latent factor models, this is computationally infeasible on large data sets. Instead, @@ -41,22 +42,23 @@ overall mean Y.. The predicted value of Yij is then Yi. + Y.j - Y.. -Computation is simple, conducted by our function **findYdotsMM**; -prediction is done on the output by our function **predict.ydotsMM**. +Computation is simple, with estimation conducted by our function +**findYdotsMM()**; prediction is done on the output by our function +**predict.ydotsMM()**. A novel enhancement in the package is to allow for different weights to be given to the αi and βj components in -the MM version. +the MM version. (With MLE it wouldn't matter, just changing the +variances.) We do make MLE available. Here αi and βj are assumed to have independent normal distributions -with different variances. (The error term -εij = Yij - EYij -is assumed independent of αi and βj, -with variance constant across i and j.) -We piggyback R's **lme4** package, -forming a wrapper for our application, and adding our function -**predict.ydotsMLE**, also a wrapper suited for our context. Since MLE +with different variances. (The error term εij = +Yij - EYij is assumed independent of +αi and βj, with variance constant +across i and j.) We piggyback R's **lme4** package, forming a wrapper +for our application, and adding our function **predict.ydotsMLE()** for +prediction, also an **lme4** wrapper suited for our context. Since MLE computation can be voluminous, our package offers a parallel version. Covariates are allowed for both the MM and MLE versions. @@ -75,7 +77,7 @@ Covariates are allowed for both the MM and MLE versions. > ivl$studage <- as.numeric(ivl$studage) > ivl$lectage <- as.numeric(ivl$lectage) > ivl$service <- as.numeric(ivl$service) -> # Make correct format, choose covs: +> # Make correct format (user ID, item ID, rating), choose covs: > ivl <- ivl[,c(1,2,7,3:6)] > # Create dummy variables in place of dept: > library(dummies) @@ -107,7 +109,7 @@ Covariates are allowed for both the MM and MLE versions. Let A denote the matrix of ratings, with Yij in row i, column j. Most of A is unknown, and we wish to predict the unknown values. -*Nonnegative matrix factorization* does this as follows: +Nonnegative Matrix Factorization (NMF) does this as follows: We find nonnegative matrices W and H, each of rank k, such that A is approximately equal to the product WH. Here k is a user-defined tuning From 57c41243ab2ed336744cb103d2226357bb645e19 Mon Sep 17 00:00:00 2001 From: ariel-shin Date: Wed, 10 May 2017 23:33:00 -0700 Subject: [PATCH 02/30] moved vignettes to inst directory --- inst/rectoolShiny2/.DS_Store | Bin 6148 -> 6148 bytes inst/vignettes/rectools.Rnw | 415 +++++++++++++++++++++++++++++++++++ 2 files changed, 415 insertions(+) create mode 100644 inst/vignettes/rectools.Rnw diff --git a/inst/rectoolShiny2/.DS_Store b/inst/rectoolShiny2/.DS_Store index 5008ddfcf53c02e82d7eee2e57c38e5672ef89f6..1a71cb0776ca1e4c312bf3287e5b925ef4942c2c 100755 GIT binary patch delta 77 zcmZoMXfc?eJ=s8n#gv61g&~olgdwpkxF|0tKQEnufq@Z-nSdAsKx`luU_h4Nn7Ei7 Lq Date: Wed, 10 May 2017 23:37:33 -0700 Subject: [PATCH 03/30] removing from main directory --- vignettes/rectools.Rnw | 415 ----------------------------------------- 1 file changed, 415 deletions(-) delete mode 100644 vignettes/rectools.Rnw diff --git a/vignettes/rectools.Rnw b/vignettes/rectools.Rnw deleted file mode 100644 index a0856f7..0000000 --- a/vignettes/rectools.Rnw +++ /dev/null @@ -1,415 +0,0 @@ -% \documentclass[a4paper,man,natbib]{apa6} -\documentclass[11pt]{article} - -\usepackage[english]{babel} -\usepackage[utf8x]{inputenc} -\usepackage{amsmath} -\usepackage{graphicx} -\usepackage[colorinlistoftodos]{todonotes} -\usepackage{listings} -\usepackage{float} -\usepackage[section]{placeins} -\setlength{\parindent}{0in} - -% library(knitr) -%\VignetteIndexEntry{Partools} - - -\title{Rectools} -% \shorttitle{Your APA6-Style Manuscript} -\author{Pooja Rajkumar, Norman Matloff} -% \affiliation{University of California, Davis} - -% \abstract{Recommendation engines have a number of different -% applications. From books to movies, they enable the analysis and -% prediction of consumer preferences. The prevalence of recommender -% systems in both the business and computational world has led to clear -% advances in prediction models over the past years. Current R packages -% include recosystem and recommenderlab. However, our new package, -% rectools, currently under development, extends its capabilities in -% several directions. One of the most important differences is that -% rectools allows users to incorporate covariates, such as age and gender, -% to improve predictive ability and better understand consumer behavior. -% Our software incorporates a number of different methods, such as -% non-negative matrix factorization, random effects models, and nearest -% neighbor methods. In addition to our incorporation of covariate -% capabilities, rectools also integrates several kinds of parallel -% computation.} - -\begin{document} -\SweaveOpts{concordance=TRUE} - -\maketitle - -{\Large UNDER CONSTRUCTION} - -\section{Recommendation Engines} - -A recommender system is an engine to predict the rating or preference -that a user would give an item. E-commerce websites in particular rely heavily on recommender systems to improve customer satstifaction, predict consumer trends, and develop products based on consumer tastes. -%which traditionally utilize collaborative filtering for their methods. -This package builds on previous methodology/packages by including several -interesting novel features, notably (a) inclusion of covariates and (b) -parallelized computation. - -\section{Methods} - -We have a data matrix of user ratings of some items. For concreteness, -let's use the common example of users rating movies. - -Though a nearest-neighbor approach is planned for addition to {\bf -rectools}, at present the package is primarily focused on what are -variously called {\it latent factor} (Koren, {\it et al}, 2009), {\it -baseline predictors} (Ekstrand {\it et al}, 2010) or {\it bias term} -(Koren {\it et al}, 2009) approaches: - -\subsection{Matrix Factorization Approach} - -Let $A$ denote the full ratings matrix, with one row for each user and -one column for each item. Thus the dimension of $A$ is $m \times n$, -where $m$ and $n$ are the numbers of users and items, respectively. - -Most of $A$ consists of unknown quantities, but we wish to find known -matrices whose product approximates $A$: - -\begin{equation} -A \approx P Q -\end{equation} - -where the matrices $P$ and $Q$ have dimensions $m \times k$ and $k -\times n$, respectively. The numbers of columns of $P$ and rows of $Q$, -$k$, are chosen to be much less than $m$ and $n$, to avoid overfitting. -In this manner, we formulate predictions for the unknown entries of $A$. - -The intuition is as follows. Consider the matrix $Q$. It has a column -for each item, so a row corresponds loosely to a set of ratings of all -the items. They are not true ratings, as they will be multiplied by -numbers in $P$, but again loosely speaking, we have found $k$ typical -user rating patterns that summarize user behavior. - -Similarly, the $k$ columns of $P$ summarize item ``behavior,'' i.e.\ how -items vary from one another. - -\subsection{ANOVA Models} - -These take a statistical Analysis of Variance approach. The rating by -user $i$ of item $j$, denoted $Y_{ij}$, is assumed to have the form - -\begin{equation} -Y_{ij} = \mu + \alpha_i + \beta_j + \epsilon_{ij} -\end{equation} - -Here $\alpha_i$ is the tendency of user $i$ to rate items higher or -lower than does the average user; $\beta_j$ is the tendency of item -$j$ to be rated higher or lower than the average item (i.e.\ this item's -relative popularity); and $\epsilon_{ij}$ is the combined effect of all -unknown factors, e.g.\ the user's mood at the time the rating is made. -The terms $\alpha_i$, $\beta_j$ and $\epsilon_{ij}$ are modeled as -random variables, so we have a {\it random effects model} in ANOVA -terminlogy, with mean 0. The parameter $\mu$ is a fixed but unknown -constant, equal to the theoretical average of all possible items by all -possible users. - -To do prediction, one uses the known user-item data to form estimates of -$\mu$ and the $\alpha_i$ and $\beta_j$. Using the statistical ``hat'' -notation to denote ``estimate of,'' our predicted rating of item $j$ -by user $i$ is - -\begin{equation} -\label{anovapred} -\widehat{Y}_{ij} = \widehat{\mu} + \widehat{\alpha}_i + -\widehat{\beta}_j -\end{equation} - -The issue then is a matter of deciding how to obtain the estimates on -the right-hand side of (\ref{anovapred}). - -Random effects ANOVA models have been used since the early years of -statistics. Typically the $\alpha_i$, $\beta_j$ and $\epsilon_{ij}$ are -assumed to have normal distributions (with different variances). The -quantities in (\ref{anovapred}) are then obtained by the Maximum Likelihood -Estimation (MLE) approach. - -The R package $lme4$ is widely used for random effects models, and is -one of the ingredients in {\bf rectools}. - -An alternative is to use another time-honored statistical technique, the -Method of Moments (MM) (Owen and Gao, 2015; Perry 2015), as follows: - -Define - -\begin{equation} -\label{yidot} -Y_{i.} = \frac{1}{N_i} \sum_{j=1}^{N_i} Y_{ij} -\end{equation} - -\begin{equation} -\label{ydotj} -Y_{.j} = \frac{1}{N_j} \sum_{i=1}^{N_i} Y_{ij} -\end{equation} - -\begin{equation} -\label{ydotdot} -Y_{..} = - \sum_{i=1}^{N_i} -\sum_{j=1}^{N_j} Y_{ij} -\end{equation} - -where $N_i$ is the number of items rated by user $i$, $N_j$ is the -number of users rating item $j$, and $N$ is the sum of all the $N_i$ -(same as sum of all the $N_j$). - -Due to the mean-0 nature of the various terms, (\ref{yidot}), (\ref{ydotj}) -and (\ref{ydotdot}) have expected values -$\alpha_i + \mu$, -$\beta_j + \mu$ and -$\alpha_i + \beta_j + \mu$, respectively. This leads to the natural -estimates - -\begin{equation} -\widehat{\alpha}_i = Y_{i.} - Y_{..} -\end{equation} - -\begin{equation} -\widehat{\beta}_j = Y_{.j} - Y_{..} -\end{equation} - -\begin{equation} -\widehat{\mu}_j = Y_{..} -\end{equation} - -These are used in (\ref{anovapred}) to obtain our actual predicted -ratigs. - - -\subsection{findYdots} - -findYdots allows us to use a latent factor model to predict values in our data set. Ydots makes use of the following equation: - - -Suppose we have a user Ali. Ali has seen the following movies: - - -\begin{tabular}{l l l} - -\textbf{userID} & \textbf{movieID} & \textbf{rating}\\ - -13 & 10& 2\\ 2 & 100 & 3 \\ .. & ... & ... \\ \end{tabular} \\ Our -estimated $\alpha_i$ here is the tendency for a user to rate a -particular item compared to everyone else. Thus, the $\alpha_{Ali}$ here -would be Ali's tendency to be an "easy" or "hard" rater. For example, if -the average rating for movies is a 3 and Ali is generally harsh grader -(with primarily 1s and 2s) then Ali's $\alpha_{Ali}$ would -be negative. Vice versa, Ali's alpha would be positive if Ali is a -generally easy grader. - -Our estimated $\beta_j$ is the tendency for the movie to -be rated really highly or poorly. For example, a movie that does really -well will have generally high ratings, and thus, a generally high -$\beta_j$. If a movie has generally poor ratings, then the -$\beta_j$ will be low. - -\subsection{Method of Moments} - -Method of Moments regresses the ratings against covariates such as age, -gender, or genre. We then subtract these predictions from the actual -value and apply the latent factor model. In order to use method of -moments, make sure your data set (called ratingsIn) follows the same -format as below: - - -\begin{tabular}{l l l} - -\textbf{userID} & \textbf{itemID} & \textbf{rating}\\ - -8 & 1& 2\\ -99 & 5& 5 \\ -.. & ... & ... \\ -\end{tabular} -\\ -\subsection{An Example} -Suppose we have a user named Karina. Suppose we have the following subsection of ratings: -\\ -\\ -\begin{tabular}{l l l} - -\textbf{userID} & \textbf{itemID} & \textbf{rating}\\ - -Karina & The Little Mermaid& 2\\ -Karina & The Avengers& 5 \\ -Karina & Captain America & 4\\ -Karina & Iron Man & ???\\ -John & Iron Man& 4\\ -Max & Iron Man & 5 \\ - -\end{tabular} -\\ -\\ - -From the looks of it, Karina seems to enjoy action/superhero movies. Given this data set, lets predict how Karina would enjoy the movie Iron Man. Thus, in this case, we want to find out our predicted: \begin{equation} Y_{Karina, Iron Man} \end{equation} -Thus, we plug into the equations above: - -\begin{equation} -\widehat{\alpha}_{Karina} = Y_{Karina.} - Y_{..} -\end{equation} - -Here, $Y_{..}$ = $\frac{2 + 5 + 4 + 4 + 5}{5}$ = 4, which is the average rating for all movies by all users. - -Here, our $Y_{Karina}$ = $\frac{2 + 5 + 4}{3}$ = 3.22 - 4 = -0.78 The 2,5, and 4 refer to Karina's ratings for The Little Mermaid, The Avengers, and Captain America respectively. After subtracting, we find that Katrina's estimated alpha is negative, thus implying that Katrina is a rather harsh grader. - - -\begin{equation} -\widehat{\beta}_{Iron Man} = Y_{.Iron Man} - Y_{..} -\end{equation} - -Next, our $\widehat{\beta}_{Iron Man}$ = $\frac{4 + 5}{2}$ = 4.5 - 4 = 0.5. The 4 and 5 refer to all the ratings for Iron man, which were given be John and Max respectively. After subtracting from overall ratings, we discover that Iron Man has a positive expected beta value. Thus, Iron Man generally received high reviews and was well liked in this data set. - - -\begin{equation} -\widehat{\mu}_{Iron Man} = Y_{..} -\end{equation} - -Lastly, $\widehat{\mu}_{Iron Man}$ is simply 4. - -In order to get our predicted value: -\begin{equation} -\widehat{Y}_{karina, Iron Man} = \widehat{\mu} + \widehat{\alpha}_{Karina} + \widehat{\beta}_{Iron Man} -\end{equation} - - -Thus, $\widehat{Y}_{karina, Iron Man}$ = 4 + (-0.78) + 0.5 = 3.72. Our predicted value for Katrina's rating for Iron Man is thus 3.72. - -\subsection{Movie Lens data (no covariates)} -Now, we will walk through a step by step example. Suppose we have a movie lens data set, which can be found here: xyz. The three columns are referred to the userID, movieID, and rating for that particular user and movie on a scale of 1-5 respectively. For the sake of convenience, we will rename the three columns like so: - -%NM ADD CODE HERE -%NM -%NM Suppose we take a snippet of the following data set: -%NM -%NM \begin{tabular}{l l l} -%NM -%NM \textbf{userID} & \textbf{movieID} & \textbf{rating}\\ -%NM -%NM 196 & 242 & 3\\ 186 & 302 & 3 \\ 22 & 377 & 1 \\ -%NM 244 & 51 & 2 \\ 166 & 346 & 1 \\ 298 & 474 & 4 \\ -%NM 115 & 265 & 2 \\ 253 & 465 & 5 \\ 305 & 451 & 3 \\ -%NM .. & .. & .. \\ -%NM -%NM \\ \end{tabular} -%NM -%NM Suppose we want to predict the rating user 196 would give movie 2. We would first make a call to findYdots. There are two seperate ways that we can predict- using either findYdotsMM or findYdotsMLE. The first uses method of moments while the later uses maximum likeihood methods. -%NM \\ Supose we set the movie lens data set under an alias, data. Looking at findYdotsMM first, we make the following call: -%NM \\ -%NM dataMeans = findYdotsMM(data) -%NM \\ -%NM This call returns the overall mean of the entire data set, the mean of all users, and the mean of all items, like so: -%NM \\ -%NM ADD PICTURE HERE -%NM \\ -%NM These means can be accessed quite easily: -%NM overall_mean = dataMeans$grandMeans -%NM user_mean = dataMeans$usrMeans -%NM item_mean = dataMeans$itmMeans -%NM \\ -%NM Note: itmMeans and usrMeans are both vectors, because they are the means of each individual user or movie. \\ -%NM Additionally, findYdotsMM has the option to simply regress the estimated mean coefficents, by signaling regressYdots like so: -%NM \\ -%NM test = findYdotsMM(data, regressYdots = TRUE) -%NM \\ -%NM This call returns all the means as well as the regressYdots attribute, which is the coefficents of the regressed means, which can be found like so:\\ -%NM coef_means - test$regressYdots \\ -%NM Which looks like so: \\ -%NM -%NM ADD PICTURE HERE -%NM -%NM \\ -%NM Additionally, findYdots has a parallelization feature for computational speed up and large data sets. \\ -%NM vals = findYdots(data, regressYdots = FALSE, cls)\\ -%NM \\ -%NM A helpful feature for users who wish to use one data set to predict another, there is a predict function for ydotsMM. Predict requires a ydotsObj and a testSet. Suppose in this case, we want to predict on the movielens data set, from our own data set: -%NM \\ -%NM testSet = data -%NM testSet$rating = 0 -%NM val = predict(dataMeans,testSet) -%NM \\ -%NM val here would simply be an array of predicted values for each user and each movie in the movielens data set, like so: -%NM \\ -%NM ENTER A PICTURE HERE -%NM \\ -%NM Suppose we wanted to predict using maximum likeihood methods:\\ -%NM res_mle = findYdotsMLE(data) \\ -%NM The means are fitted to the lme4 model and returned as a ydotsMLE type: -%NM \\ ENTER A PICTURE HERE \\ -%NM Parallelization is also a feature for lme4: \\ -%NM res_mle = findYdotsMLE(data, cls) \\ -%NM Predict is also enabled for findYdotsMLE. \\ -%NM pred_val = predict(res_mle, testSet)\\ -%NM pred_val contains a -%NM -%NM \subsection{How findYdots is called} -%NM -%NM \begin{lstlisting} -%NM findYdotsMM <- function(ratingsIn,regressYdots=FALSE,cls=NULL) -%NM findYdotsMLE <- function(ratingsIn,cls=NULL) -%NM -%NM \end{lstlisting} -%NM Both functions will return different classes. -%NM \\ -%NM findYdotsMM will return a class called ydotsMM, which contains a list of the overall mean ($Y_..$), the user means ($Y_i.$) the item means ($Y_.j$) and, if enabled, the regressed means (refer to regressYdots below). -%NM \\ -%NM findYdotsMLE will return a class called ydotsMLEpar, which contains $Y_..$, $Y_i.$, $Y_.j$. -%NM -%NM \subsection{regressYdots} -%NM -%NM If regressYdots is true, apply lm to the estimated latent factors and -%NM their product, enabling rating prediction from the resulting linear -%NM function of the factors. This is currently implemented as if there are -%NM no covariates. -%NM -%NM \subsection{Prediction} -%NM Additionally, findYdots has a predict extension on R. Thus, you can call predict.ydotsMM to call method of moments on your data set like so: -%NM predict.ydotsMM(ydotsObj,testSet), where ydotsObj is an object of type findYdots, and testSet is in the same form of ratingsIn. Similiarly, you can call predict.ydotsMLE(ydotsObj,testSet) in the same fashion. Both will return a set of predicted values. -%NM -%NM \subsection{cls} -%NM -%NM Our package is parallelized in all capacities. Thus, if a user would like to improve their computation time, a user can send in nodes to the function itself, and findydots will parallelize as needed. -%NM -%NM \subsection{xval} -%NM Additionally, you can use xvalMLE or xvalMM to split your input data into training and test sets and then determine the accuracy of your predicted values by training part of your data set to predict known values in your test set. -%NM -%NM \subsection{xvalMLE} -%NM xvalMLE can be called: -%NM \\ -%NM xvalMLE(ratingsIn, trainprop = 0.5, 'exact') -%NM \\ -%NM When using xvalMLE, you will need a data set in the same form of ratingsIn (as discussed above), with the option of adding a training proportion, an accuracy measure, and nodes to parallelize. The training proportion will be the proportion of the data set you would like to train in order to predict your test set. Your accuracy measure can be one of three different possibilties: -%NM \\ -%NM 1. exact- The accuracy is measured based on the average of all the correct predicted values compared to the test set. \\ -%NM 2. mad- The accuracy is measured based on the average of the absolute value of differences between the predicted value and actual value. \\ -%NM 3. rms - The accuracy is measured by first subtracting the predicted and actual values, squaring it, taking the mean, and then taking the square root of the mean.\\ -%NM \\ -%NM xvalMLE returns a class called result of type xvalb, which contains the number of rows in the data set, the training proportion, the accuracy measure, the number of predicted NA values, and the actual calculated accuracy. -%NM -%NM \subsection{xvalMM} -%NM xValMM is similiar to xvalMLE, except it adds a regressYdots option as discussed above. It returns a class called result of type xvalb, which contains the number of rows in your data set, the training proportion, the accuracy measure, the number of predicted NA values, and the resulting accuracy. -%NM -%NM \subsection{xvalReco} -%NM xvalReco is an option for users who want to use a parallelized version of recosystem. xvalReco requires ratingsIn, but has the option to add training proportions, accuracy measures, clusters, and a rank. xvalReco not only calls recosystem on the data set, but also seperates the data set into training and test sets and calculated the accuracy of the cross validation. xvalReco also extends the predict capability, given a recoObj and a testset. A recoObj would be the object returned from recosystem. Thus, xvalReco returns the number of rows in ratingsIn, the training proportion, the accuracy measure, the number of predicted NA values, and the accuracy.xvalReco has additional helpful functions, such as: -%NM \\ -%NM 1. getTrainSet- returns a training set given a data set and an option to add a training proportion \\ -%NM \\ 2. getTestSet- returns a testset given a data set and the training set -%NM \\ -%NM 3. trainReco- runs recosystem on a training set, with additional options to add a training proportion and a rank, and returns a class res of type Reco which contains all the predicted user and movie values. \\ -%NM \\ -%NM 4. predict- the object returned by recosystem contains two vectors, P and Q, which contain the predicted values for each user and movie respectively. Thus, in order to get the actual predicted value for a particular user and movie, predict multiplies the p and q matricies. Once these are multiplied, we will have the correct predicted values for each $Y_ij$. -%NM \\ -%NM xvalReco calls trainset, then testset, then train, and then predict. Once the predicted values are returned, xvalReco calculates the accuracy. -\section{References} - -(To be filled in.) - -Koren, 2009. Art Owen and Kaitlyn Gao, 2015. Patrick Perry, 2015. - - -\end{document} - From ae4b8c474f556808dc025a251e29ce1ea847c166 Mon Sep 17 00:00:00 2001 From: ariel-shin Date: Thu, 11 May 2017 12:03:08 -0700 Subject: [PATCH 04/30] moved findAlphBet.R.iter to inst --- {R => inst}/findAlphBet.R.iter | 0 1 file changed, 0 insertions(+), 0 deletions(-) rename {R => inst}/findAlphBet.R.iter (100%) diff --git a/R/findAlphBet.R.iter b/inst/findAlphBet.R.iter similarity index 100% rename from R/findAlphBet.R.iter rename to inst/findAlphBet.R.iter From 8ba709d2e2c8592a0d1dceb56e3b849d6bdb89ea Mon Sep 17 00:00:00 2001 From: vchakrab Date: Thu, 11 May 2017 22:55:09 -0700 Subject: [PATCH 05/30] added hybrid.R + other changes formUserDatum temporarily uncommented for testing --- .DS_Store | Bin 6148 -> 6148 bytes R/cosine.R | 1 + R/findUsrItmData.R | 14 ++--- R/hybrid.R | 126 +++++++++++++++++++++++++++++++++++++++++++++ inst/.DS_Store | Bin 6148 -> 6148 bytes 5 files changed, 134 insertions(+), 7 deletions(-) create mode 100644 R/hybrid.R diff --git a/.DS_Store b/.DS_Store index 3e280092a69b868e077e9dc20033d01d04d58e4b..bd63b964291db640c7436893bca7f8c3defcc706 100644 GIT binary patch delta 50 zcmZoMXfc=&G4YhKBr^j80}DeYLmophLrGb1QC?1dUOG^Wabl15#+!QV6C1cTvvd6A F2LNvH4+j7M delta 45 zcmZoMXfc=&F)_-JnSp_M;u39PPKGjuOonuZJcd+;5+F=vDBgI$mVGlj$6tN`6$lJV diff --git a/R/cosine.R b/R/cosine.R index c918b37..aeb4e0b 100644 --- a/R/cosine.R +++ b/R/cosine.R @@ -35,6 +35,7 @@ predict.usrData <- function(origData,newData,newItem, # (NA,NA); defined for use by sapply() below checkNewItem <- function(oneUsr) { tmp <- match(oneUsr$itms, newItem) + #message(oneUsr) if (all(is.na(tmp))) { c(NA,NA) } diff --git a/R/findUsrItmData.R b/R/findUsrItmData.R index d5402eb..11911d5 100644 --- a/R/findUsrItmData.R +++ b/R/findUsrItmData.R @@ -78,13 +78,13 @@ formUserData <- function(ratingsIn,usrCovs=NULL,itmCats=NULL,fileOut='') { retval } -### # construct a new object of class 'usrDatum' -### -### formUserDatum <- function(itms,ratings,userID=NULL) { -### obj <- list(itms = itms, ratings=ratings,userID=userID) -### class(obj) <- 'usrDatum' -### obj -### } +# construct a new object of class 'usrDatum' + +formUserDatum <- function(itms,ratings,userID=NULL) { + obj <- list(itms = itms, ratings=ratings,userID=userID) + class(obj) <- 'usrDatum' + obj +} # utility: find input row for a given user, item findInputRow <- function(ratingsIn,usrID,itmID) { diff --git a/R/hybrid.R b/R/hybrid.R new file mode 100644 index 0000000..819105a --- /dev/null +++ b/R/hybrid.R @@ -0,0 +1,126 @@ +#predict using the k nearest neighbour from the "basis" of the user-space(origData) +# covariates (e.g. age, gender) and item type preferences (e.g. +# preferred movie genres) are allowed + +########################## predict() function ########################## + +# predict newData argument from origData argument using the hybrid kNN-Matrix Product Method + +# arguments: + +# origData: training set, object of class 'usrData' (see file +# findUsrItmData.R) +# newData: data point (just one for now) to be predicted, object of +# class 'usrDatum' +# newItem: ID of the item rating to be predicted for the user in +# newData +# wtcovs: weight to put on covariates; NULL if no covs +# wtcats: weight to put on item categories; NULL if no cats +# k: a list of the number of nearest neigbhors to be considered while predicting +# r: rank for factorisation + +# value: + +# predicted ratings for newData + +hybdpredict.usrData <- function(datain,newData,newItem, + k,r, wtcovs=NULL,wtcats=NULL) { + + origData <- formUserData(datain) + #getting the "basis" matix of the "user space" + if (r >= min(dim(origData))) { + stop('no data on this item') + } + else{ + appMatrix = Reco() + train_set = data_memory(unlist(ratingsIn[1]), unlist(ratingIn[2]), + unlist(ratingsIn[3]), index1 = 1) + appMatrix$train(train_set) + #get factors after NMF -- need to use $tune() fix later + factors <- appMatrix$output(out_memory(), out_memory()) + #P- is the "user space basis"-- see RecoSys for P and Q + P <- factors$P + usrs <- rep(1:dim(P)[1], each = dim(P)[2]) + itms <- rep.int(1:dim(P)[2], dim(P)[1]) + #get ratings from "user space basis" + rtngs <- as.vector(t(P)) + #form basis (is of type usrData) + basis <- formUserData(data.frame(usrs, itms, rtngs)) + } + + # action of checkNewItem(): here oneUsr is one representative-user record in + # the "basis" of the use-space; + #the function will check which item in the record, if any, + # equals newItem, and will return that value and the corresponding + # rating; defined for use by sapply() below + checkNewItem <- function(oneUsr) { + #itms <- as.list(as.data.frame(t(oneUsr$itms))) + tmp <- match(oneUsr$itms[[1]], newItem) + if (all(is.na(tmp))) { + stop('no data on this item') + } + whichOne <- which(!is.na(tmp)) + c(whichOne,oneUsr$ratings[whichOne]) + } + found <- as.matrix(sapply(basis,checkNewItem)) + # found is of dimensions 2 x number of users; + # found[1,i] = j means origData[[i]]$itms[j] = newItem; + # found[1,i] = NA means newItem wasn't rated by user i;ß + # found[2,i] will be the rating in the non-NA case + # we need to get rid of the NA users + whoHasIt <- which(!is.na(found[1,])) + # whoHasIt[i] is the index, i.e. user ID, of the i-th user who has + # rated newData + if (is.null(whoHasIt) | length(whoHasIt)==0) return(NA) # no one rated this item + basis <- basis[whoHasIt] + # now basis only has the relevant users, the ones who have rated + # newItem, so select only those columns of the found matrix + found <- found[,whoHasIt,drop=FALSE] + #found <- found[!is.na(found)] + # find the distance from newData to one user y of origData; defined for + # use in sapply() below + onecos <- function(y) cosDist(newData,y,wtcovs,wtcats) + cosines <- sapply(basis,onecos) + # the vector cosines contains the distances from newData to all the + # original data points; + + #action of findKnghbourRtng(): predict rating based on each k[i] neighbours + #x = k[i] + #if x > neighbours present in the dataset, then the maximum number of meighbours is used + findKnghbourRtng <- function(x){ + #x can be atmost the number of neighbours in the dataset + x <- min(x, length(cosines)) + #klarge is a vector containing the indices of the x closest neighbours + klarge <- order(cosines,decreasing=TRUE)[1:x] + mean(as.numeric(found[2, klarge])) + } + sapply(k, findKnghbourRtng) +} + +# cosDist() find cosine distance between x and y, elements of an object +# of 'usrData' class; only items rated in both x and y are used; if none +# exist, then return NaN +cosDist <- function(x,y,wtcovs,wtcats) { + # rated items in common + commItms <- intersect(x$itms[[1]],y$itms[[1]]) + if (is.null(commItms)| length(commItms)==0) return(NaN) + # where are they in x and y? + xwhere <- which(!is.na(match(x$itms[[1]],commItms))) + ywhere <- which(!is.na(match(y$itms[[1]],commItms))) + xrats <- x$ratings[xwhere,] + yrats <- y$ratings[ywhere,] + cosTot <- xrats[[1]] %*% yrats[[1]] + xl2 <- sum(xrats^2) + yl2 <- sum(yrats^2) + if (!is.null(wtcovs)) { + cosTot <- cosTot + wtcovs * x$cvrs %*% y$cvrs + xl2 <- xl2 + sum((wtcovs*x$cvrs)^2) + yl2 <- yl2 + sum((wtcovs*y$cvrs)^2) + } + if (!is.null(wtcats)) { + cosTot <- cosTot + wtcats * x$cats %*% t(y$cats) + xl2 <- xl2 + sum((wtcovs*x$cats)^2) + yl2 <- yl2 + sum((wtcovs*y$cats)^2) + } + cosTot / sqrt(xl2 * yl2) +} diff --git a/inst/.DS_Store b/inst/.DS_Store index bb989f701bbbd3a1e10191f63c455fd74981aa74..c1e0b87a97c774d84bf1de3d2247fe04908720e4 100644 GIT binary patch delta 53 zcmZoMXfc@J&&a$nU^gQp^JINyxy@QkMvTIo3}pXY&nY delta 38 ucmZoMXfc@J&&aefU^nAr0}+ Date: Mon, 22 May 2017 23:02:39 -0700 Subject: [PATCH 06/30] Added buildMatrix and trainNM() --- R/findYdotsMM.R | 52 +++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 52 insertions(+) diff --git a/R/findYdotsMM.R b/R/findYdotsMM.R index abe38f6..38f4f30 100644 --- a/R/findYdotsMM.R +++ b/R/findYdotsMM.R @@ -88,4 +88,56 @@ predict.ydotsMM = function(ydotsObj,testSet) { testSet$pred = testSet$pred + predict(ydotsObj$regObj,testSet[,-(1:2)]) testSet$pred } +# buildMatrix(): A function that changes ratingsIn into a mostly empty matrix upon which we fill +buildMatrix <- function(ratingsIn,NAval=NA){ + # deal with possible factors + dmax <- function(d) { + if (is.factor(d)) return(length(levels(d))) + max(d) + } + users = ratingsIn[,1] + movies = ratingsIn[,2] + rating = ratingsIn[,3] + newMatrix = matrix(NA, + nrow = dmax(users), ncol = dmax(movies)) + for(rows in 1:nrow(ratingsIn)){ + newMatrix[ratingsIn[rows,1],ratingsIn[rows,2]] = ratingsIn[rows,3] + } + return (newMatrix) +} + +# Wrapper for recosystem +# TrainNM will check if recosystem is called, if so, call recosystem function +# If not, build ratingsIn into a mostly matrix which would be users x items dimensions +# Approximate the unknown ratings using findYdotsMM +# Substitute those values and call NMF() +trainNM <- function(ratingsIn, trainprop = 0.5,cls = NULL, + rnk = 10, recosystem = FALSE,regressYdots=FALSE){ + require(NMF) + library(NMF) + if(recosystem == TRUE){ + source((paste(system.file(package = 'rectools'), + 'recosys/xvalRecoParallel.R', sep = ""))) + res <- xvalReco(ratingsIn,trainprop,cls,rnk) + } + else { + fullMatrix <- buildMatrix(ratingsIn) # Matrix A (Step 1) + fullMatrix[which(fullMatrix == 0)] = NA + approxMatrix <- findYdotsMM(ratingsIn) # Matrix V (Step 2) + + naMatrix <- as.data.frame(which(is.na(fullMatrix) == TRUE, arr.ind = TRUE)) + naMatrix$ratings <- NA + + + preds <- predict.ydotsMM(approxMatrix,naMatrix) # Step 3 + + + + fullMatrix[which(is.na(fullMatrix))] <- preds + fullMatrix[fullMatrix < 0] <- 0 + require(NMF) + res <- nmf(as.matrix(fullMatrix),10) # Step 4 + } + res +} From 65c951b938bce3aeba2fcd22395d6e53f5246af4 Mon Sep 17 00:00:00 2001 From: Pooja-Rajkumar Date: Mon, 22 May 2017 23:03:47 -0700 Subject: [PATCH 07/30] removed recosystem call --- NAMESPACE | 1 - 1 file changed, 1 deletion(-) diff --git a/NAMESPACE b/NAMESPACE index c85ea04..48c4703 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -2,7 +2,6 @@ S3method(predict,usrData) S3method(predict,RecoS3) -importFrom("recosystem", Reco, data_memory, out_file) export(findYdotsMM,predict.ydotsMM,plot.ydotsMM,trainMM, From 9934af5786323283d39e3affa9290d2877f25b5b Mon Sep 17 00:00:00 2001 From: Pooja-Rajkumar Date: Mon, 22 May 2017 23:04:36 -0700 Subject: [PATCH 08/30] Updated to suggest recosystem --- DESCRIPTION | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/DESCRIPTION b/DESCRIPTION index fe5d837..403bc0c 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -11,7 +11,8 @@ Date: 12/12/2016 Title: Advanced Recommender System Description: Toolbox for recommender systems, including novel methods, plotting, parallel computation, and so on. Requires Recosystem of at least 0.4 -Depends: partools, lme4, recosystem, shiny +Depends: partools, lme4, shiny +Suggests: recosystem LazyLoad: no ByteCompile: yes NeedsCompilation: no From 056ce0b0c1b553fe6223146cfe1f7fafb0cdca5d Mon Sep 17 00:00:00 2001 From: matloff Date: Wed, 24 May 2017 12:59:19 -0700 Subject: [PATCH 09/30] multiplic model now moved to experimental --- {R => inst}/findAlphBet.R | 0 1 file changed, 0 insertions(+), 0 deletions(-) rename {R => inst}/findAlphBet.R (100%) diff --git a/R/findAlphBet.R b/inst/findAlphBet.R similarity index 100% rename from R/findAlphBet.R rename to inst/findAlphBet.R From f9850525994cfa6516a331b0d67158c89465b49c Mon Sep 17 00:00:00 2001 From: matloff Date: Wed, 24 May 2017 13:00:51 -0700 Subject: [PATCH 10/30] doc for multiplic --- inst/findAlphBet.R | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/inst/findAlphBet.R b/inst/findAlphBet.R index 3175b37..744c74f 100644 --- a/inst/findAlphBet.R +++ b/inst/findAlphBet.R @@ -2,9 +2,10 @@ # multiplicative random effects model: the rating Y, 0 or 1, is # generated by # -# EY_ij = q alpha_i beta_j +# EY_ij = b0 = b1 alpha_i beta_j -# for unobservable alpha, beta, and an unknown constant q +# for unobservable alpha, beta having means 1.0, +# and unknown constants b0, b1 # this is a very simple model, whose virtue is simplicity and quick # computation From 1854c635a6fe3d1a490d67c0b1f91683dabce454 Mon Sep 17 00:00:00 2001 From: matloff Date: Wed, 24 May 2017 13:17:48 -0700 Subject: [PATCH 11/30] finished moving multiplic --- NAMESPACE | 1 - {R => inst}/xvalMultiplic.R | 0 2 files changed, 1 deletion(-) rename {R => inst}/xvalMultiplic.R (100%) diff --git a/NAMESPACE b/NAMESPACE index 48c4703..1da796f 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -11,7 +11,6 @@ export(findYdotsMM,predict.ydotsMM,plot.ydotsMM,trainMM, plot.ydotsMM,plot.xvalb, getTrainSet,getTestSet,trainReco,xvalReco, focusGrp, - trainMultiplic,predict.MMmultiplic,xvalMultiplic, ratingness,covratingness) diff --git a/R/xvalMultiplic.R b/inst/xvalMultiplic.R similarity index 100% rename from R/xvalMultiplic.R rename to inst/xvalMultiplic.R From 60ce02cd0b32842549a60144d4c404522c73ce7e Mon Sep 17 00:00:00 2001 From: Norm Matloff Date: Wed, 24 May 2017 23:00:32 -0700 Subject: [PATCH 12/30] re-doing covariates --- R/newMM.R | 101 ++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 101 insertions(+) create mode 100644 R/newMM.R diff --git a/R/newMM.R b/R/newMM.R new file mode 100644 index 0000000..707675e --- /dev/null +++ b/R/newMM.R @@ -0,0 +1,101 @@ + +# EXPERIMENTAL + +# modified version of findYdotsMM() etc. +# +# in the original model +# +# Y_ij = mu + alpha_i + beta_j + eps +# +# break down alpha_i: +# +# alpha_i = sum_k phi_i x_ik + gamma_i + +# that means + +# E(Y_ij | X_ik, k = 1,...,p) = sum_k phi_i x_ik + +# so the phi can be estimated by lm() without a constant term + +# arguments: + +# ratingsIn: input data, with first cols (userID,itemID,rating, +# covariates); data frame, unless cls is non-null, in which +# case this argument is the quoted name of the distributed +# data frame + +# value: + +# S3 class of type "ydots", with components: + +# Y..: grand mean (0 if have covariates) +# Yi.: vector of mean ratings for each user +# Y.j: vector of mean ratings for each item +# regOjb: if have covariates, regression output, e.g. coefs + +findYdotsMM <- function(ratingsIn,regressYdots=FALSE) { + users <- ratingsIn[,1] + items <- ratingsIn[,2] + ratings <- ratingsIn[,3] + nms <- names(ratingsIn) + Y.. <- mean(ratings) + Yi. <- tapply(ratings,users,mean) # means of all ratings per user + Y.j <- tapply(ratings,items,mean) # means of all ratings per item + ydots <- list(grandMean=Y..,usrMeans=Yi.,itmMeans=Y.j) + haveCovs <- ncol(ratingsIn) > 3 + if (haveCovs) { + # regress, no constant term + frml <- as.formula(paste(nms[3],'~ .-1')) + lmout <- lm(frml,data=ratingsIn[,-(1:2)]) + # fits = lmout$fitted.values + # ratings = ratings - fits + # Y.. = 0 + Ni <- tapply(ratings,users,length) # means of all ratings per user + ydots$phi <- coef(lmout) + ydots$Ni <- coef(lmout) + } + class(ydots) = 'newMM' + invisible(ydots) +} + +# alias +trainNewMM <- newMM + +# test case +checkyd <- function() { + check <- data.frame( + userID <- c(1,3,2,1,2), + itemID <- c(1,1,3,2,3), + ratings <- 6:10) + print(check) + print(newMM(check)) + check$cv <- c(1,4,6,2,10) + print(check) + print(fnewMM(check)) +} + +# predict() method for the 'ydots' class +# +# arguments: +# +# testSet: data frame in same form as ratingsIn above except that there +# is no ratings column; thus covariates, if any, are shifted +# leftward one slot +# newMMObj: the output of newMM() +# minN: if Ni < minN and have covariates, use the latter +# +# returns vector of predicted values for testSet +predict.newMM = function(newMMObj,testSet,minN=NULL) { + haveCovs <- ncol(testSet) > 2 + # use of as.character() is to take advantage of row names, in case of + # future gaps in consecutive user IDs; not implemented yet + ts1 <- as.character(testSet[,1]) + tmp <- ifelse ( + newMMObj$Ni[ts1] >= minN, + newMMObj$usrMeans[ts1], + testSet$pred = newMMObj$usrMeans[as.character(testSet[,1])] + + newMMObj$itmMeans[as.character(testSet[,2])] - newMMObj$grandMean + testSet$pred +} + + From af24ec606d471d982ae5ffdaeb53236b379b7741 Mon Sep 17 00:00:00 2001 From: Norm Matloff Date: Fri, 26 May 2017 20:49:54 -0700 Subject: [PATCH 13/30] new covs code seems to work, but need more testing --- NAMESPACE | 4 ++-- R/newMM.R | 61 +++++++++++++++++++++++++++++++------------------------ 2 files changed, 37 insertions(+), 28 deletions(-) diff --git a/NAMESPACE b/NAMESPACE index 1da796f..86c1082 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -1,8 +1,7 @@ S3method(predict,usrData) S3method(predict,RecoS3) - - +S3method(predict,newMM) export(findYdotsMM,predict.ydotsMM,plot.ydotsMM,trainMM, findYdotsMLE,predict.ydotsMLE,trainMLE, @@ -11,6 +10,7 @@ export(findYdotsMM,predict.ydotsMM,plot.ydotsMM,trainMM, plot.ydotsMM,plot.xvalb, getTrainSet,getTestSet,trainReco,xvalReco, focusGrp, + findNewMM,trainNewMM, ratingness,covratingness) diff --git a/R/newMM.R b/R/newMM.R index 707675e..c535e50 100644 --- a/R/newMM.R +++ b/R/newMM.R @@ -33,7 +33,7 @@ # Y.j: vector of mean ratings for each item # regOjb: if have covariates, regression output, e.g. coefs -findYdotsMM <- function(ratingsIn,regressYdots=FALSE) { +findNewMM <- function(ratingsIn,regressYdots=FALSE) { users <- ratingsIn[,1] items <- ratingsIn[,2] ratings <- ratingsIn[,3] @@ -44,58 +44,67 @@ findYdotsMM <- function(ratingsIn,regressYdots=FALSE) { ydots <- list(grandMean=Y..,usrMeans=Yi.,itmMeans=Y.j) haveCovs <- ncol(ratingsIn) > 3 if (haveCovs) { - # regress, no constant term + # regress, no constant term; could do a weighted least squares, + # using the Ni, but since the latter is random too, not needed frml <- as.formula(paste(nms[3],'~ .-1')) lmout <- lm(frml,data=ratingsIn[,-(1:2)]) + ydots$lmout <- lmout # fits = lmout$fitted.values # ratings = ratings - fits # Y.. = 0 - Ni <- tapply(ratings,users,length) # means of all ratings per user - ydots$phi <- coef(lmout) - ydots$Ni <- coef(lmout) + Ni <- tapply(ratings,users,length) # number of ratings per user + ydots$Ni <- Ni } class(ydots) = 'newMM' invisible(ydots) } # alias -trainNewMM <- newMM - -# test case -checkyd <- function() { - check <- data.frame( - userID <- c(1,3,2,1,2), - itemID <- c(1,1,3,2,3), - ratings <- 6:10) - print(check) - print(newMM(check)) - check$cv <- c(1,4,6,2,10) - print(check) - print(fnewMM(check)) -} +trainNewMM <- findNewMM # predict() method for the 'ydots' class # -# arguments: +# arguments # # testSet: data frame in same form as ratingsIn above except that there # is no ratings column; thus covariates, if any, are shifted -# leftward one slot +# leftward one slot, i.e. userID, itemID, cov1, cov2... # newMMObj: the output of newMM() # minN: if Ni < minN and have covariates, use the latter # # returns vector of predicted values for testSet -predict.newMM = function(newMMObj,testSet,minN=NULL) { +predict.newMM = function(newMMObj,testSet,minN=0) { haveCovs <- ncol(testSet) > 2 # use of as.character() is to take advantage of row names, in case of # future gaps in consecutive user IDs; not implemented yet ts1 <- as.character(testSet[,1]) - tmp <- ifelse ( - newMMObj$Ni[ts1] >= minN, - newMMObj$usrMeans[ts1], - testSet$pred = newMMObj$usrMeans[as.character(testSet[,1])] + + if (!haveCovs) tmp <- newMMObj$usrMeans[ts1] else { + tmp <- ifelse ( + newMMObj$Ni[ts1] >= minN, + newMMObj$usrMeans[ts1], + predict(newMMObj$lmout,testSet[,-(1:2),drop=FALSE])) + } + testSet$pred <- tmp + newMMObj$itmMeans[as.character(testSet[,2])] - newMMObj$grandMean testSet$pred } +# test case +checkydNew <- function() { + check <- data.frame( + userID <- c(1,3,2,1,2), + itemID <- c(1,1,3,2,3), + ratings <- 6:10) + names(check) <- c('u','i','r') + print(check) + print(findNewMM(check)) + check$cv <- c(1,4,6,2,10) + names(check)[4] <- 'x' + print(check) + cout <- findNewMM(check) + print(cout) + testset <- check[1:2,-3] + testset$x <- c(5,8) + print(predict(cout,testset,2)) +} From 71ab8cd3b6abfdfefb321307bee4bfd111979641 Mon Sep 17 00:00:00 2001 From: Norman Matloff Date: Sat, 27 May 2017 15:43:34 -0700 Subject: [PATCH 14/30] almost ready --- R/newMM.R | 50 +++++++++++++++++++++++++------------------------- 1 file changed, 25 insertions(+), 25 deletions(-) diff --git a/R/newMM.R b/R/newMM.R index c535e50..5abdaa4 100644 --- a/R/newMM.R +++ b/R/newMM.R @@ -1,9 +1,8 @@ -# EXPERIMENTAL +# modified version of earlier findYdotsMM() etc., May 27, 2017, with +# new approach to use of covariates -# modified version of findYdotsMM() etc. -# -# in the original model +# in the basic model # # Y_ij = mu + alpha_i + beta_j + eps # @@ -17,12 +16,12 @@ # so the phi can be estimated by lm() without a constant term +# the same could be done for the items, not yet implemented + # arguments: -# ratingsIn: input data, with first cols (userID,itemID,rating, -# covariates); data frame, unless cls is non-null, in which -# case this argument is the quoted name of the distributed -# data frame +# ratingsIn: input data, with cols (userID,itemID,rating, +# covariates); data framee # value: @@ -49,43 +48,44 @@ findNewMM <- function(ratingsIn,regressYdots=FALSE) { frml <- as.formula(paste(nms[3],'~ .-1')) lmout <- lm(frml,data=ratingsIn[,-(1:2)]) ydots$lmout <- lmout - # fits = lmout$fitted.values - # ratings = ratings - fits - # Y.. = 0 Ni <- tapply(ratings,users,length) # number of ratings per user ydots$Ni <- Ni } - class(ydots) = 'newMM' + class(ydots) = 'ydotsMM' invisible(ydots) } # alias -trainNewMM <- findNewMM +trainMM <- findMM + +# predict() method for the 'ydotsMM' class + +# in predicting for user i, the code looks at N_i, the number of ratings +# by user i; if that number is below minN, the prediction comes from +# user i's covariate information -# predict() method for the 'ydots' class -# # arguments -# + # testSet: data frame in same form as ratingsIn above except that there # is no ratings column; thus covariates, if any, are shifted # leftward one slot, i.e. userID, itemID, cov1, cov2... -# newMMObj: the output of newMM() -# minN: if Ni < minN and have covariates, use the latter -# +# ydotsOjb: the output of findYdotsMM() +# minN: if Ni < minN and have covariates, use the latter instead of # Yi. + # returns vector of predicted values for testSet -predict.newMM = function(newMMObj,testSet,minN=0) { +predict.ydotsMM = function(ydotsOjb,testSet,minN=0) { haveCovs <- ncol(testSet) > 2 # use of as.character() is to take advantage of row names, in case of # future gaps in consecutive user IDs; not implemented yet ts1 <- as.character(testSet[,1]) - if (!haveCovs) tmp <- newMMObj$usrMeans[ts1] else { + if (!haveCovs) tmp <- ydotsOjb$usrMeans[ts1] else { tmp <- ifelse ( - newMMObj$Ni[ts1] >= minN, - newMMObj$usrMeans[ts1], - predict(newMMObj$lmout,testSet[,-(1:2),drop=FALSE])) + ydotsOjb$Ni[ts1] >= minN, + ydotsOjb$usrMeans[ts1], + predict(ydotsOjb$lmout,testSet[,-(1:2),drop=FALSE])) } testSet$pred <- tmp + - newMMObj$itmMeans[as.character(testSet[,2])] - newMMObj$grandMean + ydotsOjb$itmMeans[as.character(testSet[,2])] - ydotsOjb$grandMean testSet$pred } From 8fe40775af775f8781de7426e9c8edf51e3ba356 Mon Sep 17 00:00:00 2001 From: Norman Matloff Date: Sat, 27 May 2017 16:49:18 -0700 Subject: [PATCH 15/30] improved MM w/ covs ready --- R/findYdotsMM.R | 195 +++++++++++++++++++--------------------------- R/findYdotsMM.old | 143 ++++++++++++++++++++++++++++++++++ R/newMM.R | 110 -------------------------- 3 files changed, 224 insertions(+), 224 deletions(-) create mode 100644 R/findYdotsMM.old delete mode 100644 R/newMM.R diff --git a/R/findYdotsMM.R b/R/findYdotsMM.R index 38f4f30..5abdaa4 100644 --- a/R/findYdotsMM.R +++ b/R/findYdotsMM.R @@ -1,15 +1,27 @@ +# modified version of earlier findYdotsMM() etc., May 27, 2017, with +# new approach to use of covariates + +# in the basic model +# +# Y_ij = mu + alpha_i + beta_j + eps +# +# break down alpha_i: +# +# alpha_i = sum_k phi_i x_ik + gamma_i + +# that means + +# E(Y_ij | X_ik, k = 1,...,p) = sum_k phi_i x_ik + +# so the phi can be estimated by lm() without a constant term + +# the same could be done for the items, not yet implemented + # arguments: -# ratingsIn: input data, with first cols (userID,itemID,rating, -# covariates); data frame, unless cls is non-null, in which -# case this argument is the quoted name of the distributed -# data frame -# regressYdots; if TRUE, apply lm() to the estimated latent factors -# (and their product), enabling rating prediction from the -# resulting linear function of the factors; currently -# only implemented if have no covariates -# cls: an R 'parallel' cluster +# ratingsIn: input data, with cols (userID,itemID,rating, +# covariates); data framee # value: @@ -20,124 +32,79 @@ # Y.j: vector of mean ratings for each item # regOjb: if have covariates, regression output, e.g. coefs -findYdotsMM <- function(ratingsIn,regressYdots=FALSE,cls=NULL) { - require(partools) - if (is.null(cls)) { - users = ratingsIn[,1] - items = ratingsIn[,2] - ratings = ratingsIn[,3] - nms <- names(ratingsIn) - haveCovs = ncol(ratingsIn) > 3 - if (haveCovs) { - frml = as.formula(paste(nms[3],'~ .')) - lmout = lm(frml,data=ratingsIn[,-(1:2)]) - fits = lmout$fitted.values - ratings = ratings - fits - Y.. = 0 - } else Y.. = mean(ratings) - Yi. = tapply(ratings,users,mean) # means of all ratings per user - Y.j = tapply(ratings,items,mean) # means of all ratings per item - } else { - stop('parallel version under construction') - } - ydots = list(grandMean=Y..,usrMeans=Yi.,itmMeans=Y.j) +findNewMM <- function(ratingsIn,regressYdots=FALSE) { + users <- ratingsIn[,1] + items <- ratingsIn[,2] + ratings <- ratingsIn[,3] + nms <- names(ratingsIn) + Y.. <- mean(ratings) + Yi. <- tapply(ratings,users,mean) # means of all ratings per user + Y.j <- tapply(ratings,items,mean) # means of all ratings per item + ydots <- list(grandMean=Y..,usrMeans=Yi.,itmMeans=Y.j) + haveCovs <- ncol(ratingsIn) > 3 if (haveCovs) { - ydots$regObj = lmout + # regress, no constant term; could do a weighted least squares, + # using the Ni, but since the latter is random too, not needed + frml <- as.formula(paste(nms[3],'~ .-1')) + lmout <- lm(frml,data=ratingsIn[,-(1:2)]) + ydots$lmout <- lmout + Ni <- tapply(ratings,users,length) # number of ratings per user + ydots$Ni <- Ni } - if (regressYdots && !haveCovs) { - yi. = Yi.[ratingsIn[,1]] - y.j = Y.j[ratingsIn[,2]] - # yij = yi. * y.j - # ydots$regressYdots = coef(lm(ratings ~ yi. + y.j + yij)) - ydots$regressYdots = coef(lm(ratings ~ yi. + y.j)) - } - class(ydots) = 'ydotsMM' invisible(ydots) } # alias -trainMM <- findYdotsMM +trainMM <- findMM -# check -checkyd <- function() { - check <- - data.frame(userID = c(1,3,2,1,2),itemID = c(1,1,3,2,3),ratings=6:10) - print(check) - print(findYdotsMM(check)) - check$cv <- c(1,4,6,2,10) - print(check) - print(findYdotsMM(check)) -} +# predict() method for the 'ydotsMM' class -# predict from output from coef(lm()); not used currently -predict.lmcoef = function(coefs,testSet) { - cbind(1,as.matrix(testSet)) %*% as.vector(coefs) -} +# in predicting for user i, the code looks at N_i, the number of ratings +# by user i; if that number is below minN, the prediction comes from +# user i's covariate information + +# arguments + +# testSet: data frame in same form as ratingsIn above except that there +# is no ratings column; thus covariates, if any, are shifted +# leftward one slot, i.e. userID, itemID, cov1, cov2... +# ydotsOjb: the output of findYdotsMM() +# minN: if Ni < minN and have covariates, use the latter instead of # Yi. -# predict() method for the 'ydots' class -# -# testSet in same form as ratingsIn in findYdots(), except that there -# is no ratings column; regObj is as in the output of findYdots() -# # returns vector of predicted values for testSet -predict.ydotsMM = function(ydotsObj,testSet) { - testSet$pred = ydotsObj$usrMeans[as.character(testSet[,1])] + - ydotsObj$itmMeans[as.character(testSet[,2])] - ydotsObj$grandMean - if (!is.null(ydotsObj$regObj)) - testSet$pred = testSet$pred + predict(ydotsObj$regObj,testSet[,-(1:2)]) +predict.ydotsMM = function(ydotsOjb,testSet,minN=0) { + haveCovs <- ncol(testSet) > 2 + # use of as.character() is to take advantage of row names, in case of + # future gaps in consecutive user IDs; not implemented yet + ts1 <- as.character(testSet[,1]) + if (!haveCovs) tmp <- ydotsOjb$usrMeans[ts1] else { + tmp <- ifelse ( + ydotsOjb$Ni[ts1] >= minN, + ydotsOjb$usrMeans[ts1], + predict(ydotsOjb$lmout,testSet[,-(1:2),drop=FALSE])) + } + testSet$pred <- tmp + + ydotsOjb$itmMeans[as.character(testSet[,2])] - ydotsOjb$grandMean testSet$pred } -# buildMatrix(): A function that changes ratingsIn into a mostly empty matrix upon which we fill -buildMatrix <- function(ratingsIn,NAval=NA){ - # deal with possible factors - dmax <- function(d) { - if (is.factor(d)) return(length(levels(d))) - max(d) - } - users = ratingsIn[,1] - movies = ratingsIn[,2] - rating = ratingsIn[,3] - newMatrix = matrix(NA, - nrow = dmax(users), ncol = dmax(movies)) - for(rows in 1:nrow(ratingsIn)){ - newMatrix[ratingsIn[rows,1],ratingsIn[rows,2]] = ratingsIn[rows,3] - } - return (newMatrix) -} -# Wrapper for recosystem -# TrainNM will check if recosystem is called, if so, call recosystem function -# If not, build ratingsIn into a mostly matrix which would be users x items dimensions -# Approximate the unknown ratings using findYdotsMM -# Substitute those values and call NMF() -trainNM <- function(ratingsIn, trainprop = 0.5,cls = NULL, - rnk = 10, recosystem = FALSE,regressYdots=FALSE){ - require(NMF) - library(NMF) - if(recosystem == TRUE){ - source((paste(system.file(package = 'rectools'), - 'recosys/xvalRecoParallel.R', sep = ""))) - res <- xvalReco(ratingsIn,trainprop,cls,rnk) - } - else { - fullMatrix <- buildMatrix(ratingsIn) # Matrix A (Step 1) - fullMatrix[which(fullMatrix == 0)] = NA - approxMatrix <- findYdotsMM(ratingsIn) # Matrix V (Step 2) - - naMatrix <- as.data.frame(which(is.na(fullMatrix) == TRUE, arr.ind = TRUE)) - naMatrix$ratings <- NA - - - preds <- predict.ydotsMM(approxMatrix,naMatrix) # Step 3 - - - - fullMatrix[which(is.na(fullMatrix))] <- preds - fullMatrix[fullMatrix < 0] <- 0 - require(NMF) - res <- nmf(as.matrix(fullMatrix),10) # Step 4 - } - res +# test case +checkydNew <- function() { + check <- data.frame( + userID <- c(1,3,2,1,2), + itemID <- c(1,1,3,2,3), + ratings <- 6:10) + names(check) <- c('u','i','r') + print(check) + print(findNewMM(check)) + check$cv <- c(1,4,6,2,10) + names(check)[4] <- 'x' + print(check) + cout <- findNewMM(check) + print(cout) + testset <- check[1:2,-3] + testset$x <- c(5,8) + print(predict(cout,testset,2)) } diff --git a/R/findYdotsMM.old b/R/findYdotsMM.old new file mode 100644 index 0000000..38f4f30 --- /dev/null +++ b/R/findYdotsMM.old @@ -0,0 +1,143 @@ + +# arguments: + +# ratingsIn: input data, with first cols (userID,itemID,rating, +# covariates); data frame, unless cls is non-null, in which +# case this argument is the quoted name of the distributed +# data frame +# regressYdots; if TRUE, apply lm() to the estimated latent factors +# (and their product), enabling rating prediction from the +# resulting linear function of the factors; currently +# only implemented if have no covariates +# cls: an R 'parallel' cluster + +# value: + +# S3 class of type "ydots", with components: + +# Y..: grand mean (0 if have covariates) +# Yi.: vector of mean ratings for each user +# Y.j: vector of mean ratings for each item +# regOjb: if have covariates, regression output, e.g. coefs + +findYdotsMM <- function(ratingsIn,regressYdots=FALSE,cls=NULL) { + require(partools) + if (is.null(cls)) { + users = ratingsIn[,1] + items = ratingsIn[,2] + ratings = ratingsIn[,3] + nms <- names(ratingsIn) + haveCovs = ncol(ratingsIn) > 3 + if (haveCovs) { + frml = as.formula(paste(nms[3],'~ .')) + lmout = lm(frml,data=ratingsIn[,-(1:2)]) + fits = lmout$fitted.values + ratings = ratings - fits + Y.. = 0 + } else Y.. = mean(ratings) + Yi. = tapply(ratings,users,mean) # means of all ratings per user + Y.j = tapply(ratings,items,mean) # means of all ratings per item + } else { + stop('parallel version under construction') + } + ydots = list(grandMean=Y..,usrMeans=Yi.,itmMeans=Y.j) + if (haveCovs) { + ydots$regObj = lmout + } + if (regressYdots && !haveCovs) { + yi. = Yi.[ratingsIn[,1]] + y.j = Y.j[ratingsIn[,2]] + # yij = yi. * y.j + # ydots$regressYdots = coef(lm(ratings ~ yi. + y.j + yij)) + ydots$regressYdots = coef(lm(ratings ~ yi. + y.j)) + } + + class(ydots) = 'ydotsMM' + invisible(ydots) +} + +# alias +trainMM <- findYdotsMM + +# check +checkyd <- function() { + check <- + data.frame(userID = c(1,3,2,1,2),itemID = c(1,1,3,2,3),ratings=6:10) + print(check) + print(findYdotsMM(check)) + check$cv <- c(1,4,6,2,10) + print(check) + print(findYdotsMM(check)) +} + +# predict from output from coef(lm()); not used currently +predict.lmcoef = function(coefs,testSet) { + cbind(1,as.matrix(testSet)) %*% as.vector(coefs) +} + +# predict() method for the 'ydots' class +# +# testSet in same form as ratingsIn in findYdots(), except that there +# is no ratings column; regObj is as in the output of findYdots() +# +# returns vector of predicted values for testSet +predict.ydotsMM = function(ydotsObj,testSet) { + testSet$pred = ydotsObj$usrMeans[as.character(testSet[,1])] + + ydotsObj$itmMeans[as.character(testSet[,2])] - ydotsObj$grandMean + if (!is.null(ydotsObj$regObj)) + testSet$pred = testSet$pred + predict(ydotsObj$regObj,testSet[,-(1:2)]) + testSet$pred +} +# buildMatrix(): A function that changes ratingsIn into a mostly empty matrix upon which we fill +buildMatrix <- function(ratingsIn,NAval=NA){ + # deal with possible factors + dmax <- function(d) { + if (is.factor(d)) return(length(levels(d))) + max(d) + } + users = ratingsIn[,1] + movies = ratingsIn[,2] + rating = ratingsIn[,3] + newMatrix = matrix(NA, + nrow = dmax(users), ncol = dmax(movies)) + for(rows in 1:nrow(ratingsIn)){ + newMatrix[ratingsIn[rows,1],ratingsIn[rows,2]] = ratingsIn[rows,3] + } + return (newMatrix) +} + +# Wrapper for recosystem +# TrainNM will check if recosystem is called, if so, call recosystem function +# If not, build ratingsIn into a mostly matrix which would be users x items dimensions +# Approximate the unknown ratings using findYdotsMM +# Substitute those values and call NMF() +trainNM <- function(ratingsIn, trainprop = 0.5,cls = NULL, + rnk = 10, recosystem = FALSE,regressYdots=FALSE){ + require(NMF) + library(NMF) + if(recosystem == TRUE){ + source((paste(system.file(package = 'rectools'), + 'recosys/xvalRecoParallel.R', sep = ""))) + res <- xvalReco(ratingsIn,trainprop,cls,rnk) + } + else { + fullMatrix <- buildMatrix(ratingsIn) # Matrix A (Step 1) + fullMatrix[which(fullMatrix == 0)] = NA + approxMatrix <- findYdotsMM(ratingsIn) # Matrix V (Step 2) + + naMatrix <- as.data.frame(which(is.na(fullMatrix) == TRUE, arr.ind = TRUE)) + naMatrix$ratings <- NA + + + preds <- predict.ydotsMM(approxMatrix,naMatrix) # Step 3 + + + + fullMatrix[which(is.na(fullMatrix))] <- preds + fullMatrix[fullMatrix < 0] <- 0 + require(NMF) + res <- nmf(as.matrix(fullMatrix),10) # Step 4 + } + res +} + diff --git a/R/newMM.R b/R/newMM.R deleted file mode 100644 index 5abdaa4..0000000 --- a/R/newMM.R +++ /dev/null @@ -1,110 +0,0 @@ - -# modified version of earlier findYdotsMM() etc., May 27, 2017, with -# new approach to use of covariates - -# in the basic model -# -# Y_ij = mu + alpha_i + beta_j + eps -# -# break down alpha_i: -# -# alpha_i = sum_k phi_i x_ik + gamma_i - -# that means - -# E(Y_ij | X_ik, k = 1,...,p) = sum_k phi_i x_ik - -# so the phi can be estimated by lm() without a constant term - -# the same could be done for the items, not yet implemented - -# arguments: - -# ratingsIn: input data, with cols (userID,itemID,rating, -# covariates); data framee - -# value: - -# S3 class of type "ydots", with components: - -# Y..: grand mean (0 if have covariates) -# Yi.: vector of mean ratings for each user -# Y.j: vector of mean ratings for each item -# regOjb: if have covariates, regression output, e.g. coefs - -findNewMM <- function(ratingsIn,regressYdots=FALSE) { - users <- ratingsIn[,1] - items <- ratingsIn[,2] - ratings <- ratingsIn[,3] - nms <- names(ratingsIn) - Y.. <- mean(ratings) - Yi. <- tapply(ratings,users,mean) # means of all ratings per user - Y.j <- tapply(ratings,items,mean) # means of all ratings per item - ydots <- list(grandMean=Y..,usrMeans=Yi.,itmMeans=Y.j) - haveCovs <- ncol(ratingsIn) > 3 - if (haveCovs) { - # regress, no constant term; could do a weighted least squares, - # using the Ni, but since the latter is random too, not needed - frml <- as.formula(paste(nms[3],'~ .-1')) - lmout <- lm(frml,data=ratingsIn[,-(1:2)]) - ydots$lmout <- lmout - Ni <- tapply(ratings,users,length) # number of ratings per user - ydots$Ni <- Ni - } - class(ydots) = 'ydotsMM' - invisible(ydots) -} - -# alias -trainMM <- findMM - -# predict() method for the 'ydotsMM' class - -# in predicting for user i, the code looks at N_i, the number of ratings -# by user i; if that number is below minN, the prediction comes from -# user i's covariate information - -# arguments - -# testSet: data frame in same form as ratingsIn above except that there -# is no ratings column; thus covariates, if any, are shifted -# leftward one slot, i.e. userID, itemID, cov1, cov2... -# ydotsOjb: the output of findYdotsMM() -# minN: if Ni < minN and have covariates, use the latter instead of # Yi. - -# returns vector of predicted values for testSet -predict.ydotsMM = function(ydotsOjb,testSet,minN=0) { - haveCovs <- ncol(testSet) > 2 - # use of as.character() is to take advantage of row names, in case of - # future gaps in consecutive user IDs; not implemented yet - ts1 <- as.character(testSet[,1]) - if (!haveCovs) tmp <- ydotsOjb$usrMeans[ts1] else { - tmp <- ifelse ( - ydotsOjb$Ni[ts1] >= minN, - ydotsOjb$usrMeans[ts1], - predict(ydotsOjb$lmout,testSet[,-(1:2),drop=FALSE])) - } - testSet$pred <- tmp + - ydotsOjb$itmMeans[as.character(testSet[,2])] - ydotsOjb$grandMean - testSet$pred -} - -# test case -checkydNew <- function() { - check <- data.frame( - userID <- c(1,3,2,1,2), - itemID <- c(1,1,3,2,3), - ratings <- 6:10) - names(check) <- c('u','i','r') - print(check) - print(findNewMM(check)) - check$cv <- c(1,4,6,2,10) - names(check)[4] <- 'x' - print(check) - cout <- findNewMM(check) - print(cout) - testset <- check[1:2,-3] - testset$x <- c(5,8) - print(predict(cout,testset,2)) -} - From f3a7dc4219c06934735fd2f9f2fc2e625a529a2c Mon Sep 17 00:00:00 2001 From: Norman Matloff Date: Sat, 27 May 2017 16:55:11 -0700 Subject: [PATCH 16/30] updating pkg for new covs method --- NAMESPACE | 3 +-- R/findYdotsMM.R | 15 ++++++++------- 2 files changed, 9 insertions(+), 9 deletions(-) diff --git a/NAMESPACE b/NAMESPACE index 86c1082..cb2bde4 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -1,7 +1,7 @@ S3method(predict,usrData) S3method(predict,RecoS3) -S3method(predict,newMM) +S3method(predict,ydotsMM) export(findYdotsMM,predict.ydotsMM,plot.ydotsMM,trainMM, findYdotsMLE,predict.ydotsMLE,trainMLE, @@ -10,7 +10,6 @@ export(findYdotsMM,predict.ydotsMM,plot.ydotsMM,trainMM, plot.ydotsMM,plot.xvalb, getTrainSet,getTestSet,trainReco,xvalReco, focusGrp, - findNewMM,trainNewMM, ratingness,covratingness) diff --git a/R/findYdotsMM.R b/R/findYdotsMM.R index 5abdaa4..ad5d2db 100644 --- a/R/findYdotsMM.R +++ b/R/findYdotsMM.R @@ -32,7 +32,7 @@ # Y.j: vector of mean ratings for each item # regOjb: if have covariates, regression output, e.g. coefs -findNewMM <- function(ratingsIn,regressYdots=FALSE) { +findYdotsMM <- function(ratingsIn,regressYdots=FALSE) { users <- ratingsIn[,1] items <- ratingsIn[,2] ratings <- ratingsIn[,3] @@ -56,13 +56,13 @@ findNewMM <- function(ratingsIn,regressYdots=FALSE) { } # alias -trainMM <- findMM +trainMM <- findYdotsMM # predict() method for the 'ydotsMM' class # in predicting for user i, the code looks at N_i, the number of ratings # by user i; if that number is below minN, the prediction comes from -# user i's covariate information +# user i's covariate information instead of from Yi. # arguments @@ -82,7 +82,8 @@ predict.ydotsMM = function(ydotsOjb,testSet,minN=0) { tmp <- ifelse ( ydotsOjb$Ni[ts1] >= minN, ydotsOjb$usrMeans[ts1], - predict(ydotsOjb$lmout,testSet[,-(1:2),drop=FALSE])) + predict(ydotsOjb$lmout,testSet[,-(1:2),drop=FALSE]) + ) } testSet$pred <- tmp + ydotsOjb$itmMeans[as.character(testSet[,2])] - ydotsOjb$grandMean @@ -90,18 +91,18 @@ predict.ydotsMM = function(ydotsOjb,testSet,minN=0) { } # test case -checkydNew <- function() { +checkyd <- function() { check <- data.frame( userID <- c(1,3,2,1,2), itemID <- c(1,1,3,2,3), ratings <- 6:10) names(check) <- c('u','i','r') print(check) - print(findNewMM(check)) + print(findYdotsMM(check)) check$cv <- c(1,4,6,2,10) names(check)[4] <- 'x' print(check) - cout <- findNewMM(check) + cout <- findYdotsMM(check) print(cout) testset <- check[1:2,-3] testset$x <- c(5,8) From 113148915bcb6a28887c2e3fb32db8a6000ba078 Mon Sep 17 00:00:00 2001 From: Norman Matloff Date: Sat, 27 May 2017 17:18:53 -0700 Subject: [PATCH 17/30] updating for new covs method --- R/xvalMM.R | 11 +++++------ man/xvalReco.Rd | 6 +++++- 2 files changed, 10 insertions(+), 7 deletions(-) diff --git a/R/xvalMM.R b/R/xvalMM.R index 6bb0beb..8e83385 100644 --- a/R/xvalMM.R +++ b/R/xvalMM.R @@ -16,16 +16,15 @@ # (and their product), enabling rating prediction from the # resulting linear function of the factors; currently # only implemented if have no covariates -# cls: if non-null, +# minN: see findYdotsMM.R # value: # accuracy value xvalMM <- function(ratingsIn, trainprop=0.5, - regressYdots=FALSE,cls=NULL){ - if(!is.null(cls)) stop('parallel version under construction') - if(is.null(cls)) ratIn = ratingsIn else ratIn = get(ratingsIn) + regressYdots=FALSE,minN=0){ + ratIn = ratingsIn # split into random training and validation sets nrowRatIn = nrow(ratIn) rowNum = floor(trainprop * nrowRatIn) @@ -35,13 +34,13 @@ xvalMM <- function(ratingsIn, trainprop=0.5, trainItems = trainingSet[,2] trainUsers = trainingSet[,1] # get means - means = findYdotsMM(trainingSet,regressYdots,cls) + means = findYdotsMM(trainingSet,regressYdots) # Y.. = means$grandMean # Yi. = means$usrMeans # Y.j = means$itmMeans testIdxs = setdiff(1:nrowRatIn,trainIdxs) testA = ratIn[testIdxs,] - testA$pred = predict(means,testA[,-3]) # predict.ydotsMM + testA$pred = predict(means,testA[,-3],minN) # predict.ydotsMM # need to integrate the following into predict.ydotsMM if (regressYdots) { yi. = means$usrMeans[testA[,1]] diff --git a/man/xvalReco.Rd b/man/xvalReco.Rd index b09d4f3..e6409e8 100644 --- a/man/xvalReco.Rd +++ b/man/xvalReco.Rd @@ -16,7 +16,7 @@ computation capability in some cases. } \usage{ - xvalMM(ratingsIn, trainprop = 0.5, regressYdots = FALSE, cls = NULL) + xvalMM(ratingsIn, trainprop = 0.5, regressYdots = FALSE, minN = 0) xvalMLE(ratingsIn, trainprop = 0.5, cls = NULL) xvalReco(ratingsIn, trainprop = 0.5,cls = NULL,rnk = 10) xvalCos(ratingsIn, k, usrCovs = NULL, itmCats = NULL, wtcovs = NULL, @@ -30,6 +30,10 @@ computation capability in some cases. \arguments{ \item{ratingsIn}{Input data frame. Within-row format is UserID, ItemID, rating} + \item{minN}{Applies to situations in which covariates are + present. In predicting for user i, then either Yi. or the + regression-based prediction will be used, depending on whether + N_i >= minN.} \item{trainprop}{The fraction of ratingsIn we want to use for our training set} \item{cls}{R \code{paralell} cluster.} From 2648ea40a891ce2cf943195e26c7110d198b1429 Mon Sep 17 00:00:00 2001 From: Norm Matloff Date: Sat, 27 May 2017 23:43:22 -0700 Subject: [PATCH 18/30] easy prep for 2 data sets --- DESCRIPTION | 6 ++--- NAMESPACE | 3 ++- R/Utils.R | 65 +++++++++++++++++++++++++++++++++++++++++++++++++++++ 3 files changed, 70 insertions(+), 4 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 403bc0c..fb5fbf8 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -7,11 +7,11 @@ Author@R: c(person("Ariel" "Shin", role-c("aut"), email="arishin@ucdavis.edu") ) Maintainer: Norm Matloff -Date: 12/12/2016 +Date: 5/27/2017 Title: Advanced Recommender System Description: Toolbox for recommender systems, including novel methods, - plotting, parallel computation, and so on. Requires Recosystem of at least 0.4 -Depends: partools, lme4, shiny + plotting, parallel computation, and so on. +Depends: partools, lme4, shiny, dummies Suggests: recosystem LazyLoad: no ByteCompile: yes diff --git a/NAMESPACE b/NAMESPACE index cb2bde4..498428d 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -10,6 +10,7 @@ export(findYdotsMM,predict.ydotsMM,plot.ydotsMM,trainMM, plot.ydotsMM,plot.xvalb, getTrainSet,getTestSet,trainReco,xvalReco, focusGrp, - ratingness,covratingness) + ratingness,covratingness, + gassign,getInstEval,getML) diff --git a/R/Utils.R b/R/Utils.R index 108485e..7f89322 100644 --- a/R/Utils.R +++ b/R/Utils.R @@ -17,3 +17,68 @@ covratingness <- function(ratingsIn,splitCol) { tmp[ratingsIn[,splitCol]] } +# get the instructor evaluation data, and set it up, including +# covariates + +getInstEval <- function() { + data(InstEval) + ivl <- InstEval + # convert from factors + ivl$s <- as.numeric(ivl$s) + ivl$d <- as.numeric(ivl$d) + ivl$studage <- as.numeric(ivl$studage) + ivl$lectage <- as.numeric(ivl$lectage) + ivl$service <- as.numeric(ivl$service) + # make correct format, choose + ivl <- ivl[,c(1,2,7,3:6)] + # create dummy variables in place of dept + library(dummies) + dms <- dummy(ivl$dept) + dms <- as.data.frame(dms) + dms$dept2 <- NULL + ivl$dept <- NULL + ivl <- cbind(ivl,dms) + gassign('ivl','ivl') +} + +# get the MovieLense 100K evaluation data, and set it up, including +# covariates + +# if needDownload, then download will be done; datadir is the directory +# containing the data, assumed by default to be in the current working +# directory + +getML <- function(needDownload=FALSE,datadir='./ml-100k') { + if (needDownload) { + # 5 Mb + download.file( + 'http://files.grouplens.org/datasets/movielens/ml-100k.zip', + 'ml-100k.zip') + unzip('ml-100k.zip') + } + currdir <- getwd() # leave a trail of bread crumbs + setwd(datadir) + ud <- read.table('u.data',header=F,sep='\t') + uu <- read.table('u.user',header=F,sep='|') + ui <- read.table('u.item',header=F,sep='|') + setwd(currdir) # follow the trail back + ud <- ud[,-4] # remove timestamp, leaving user, item, rating + uu <- uu[,1:3] # user, age, gender + ui <- ui[,c(1,6:24)] # item num, genres + names(ud) <- c('user','item','rating') + names(uu) <- c('user','age','gender') + names(ui)[1] <- 'item' + names(ui)[-1] <- gsub('V','GNR',names(ui)[-1]) # GNR = genre + uu$gender <- as.integer(uu$gender == 'M') + uduu <- merge(ud,uu) + uduuui <- merge(uduu,ui) + gassign('uduu','uduu') + gassign('uduuui','uduuui') +} + +# not-so-stealth global assignment +gassign <- function(gname,rhsname) { + cmd <- paste(gname,'<<-',rhsname) + eval(parse(text=cmd),envir=parent.frame()) +} + From c3c989bfce458600f7c4afe9f5cfc1db16935d61 Mon Sep 17 00:00:00 2001 From: Norm Matloff Date: Sun, 28 May 2017 10:03:01 -0700 Subject: [PATCH 19/30] correct for gaps in IDs --- R/findYdotsMM.R | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/R/findYdotsMM.R b/R/findYdotsMM.R index ad5d2db..d275a41 100644 --- a/R/findYdotsMM.R +++ b/R/findYdotsMM.R @@ -35,6 +35,8 @@ findYdotsMM <- function(ratingsIn,regressYdots=FALSE) { users <- ratingsIn[,1] items <- ratingsIn[,2] + # user and item IDs may not be consecutive; even if they are + # consecutive in the original ratings <- ratingsIn[,3] nms <- names(ratingsIn) Y.. <- mean(ratings) @@ -70,13 +72,12 @@ trainMM <- findYdotsMM # is no ratings column; thus covariates, if any, are shifted # leftward one slot, i.e. userID, itemID, cov1, cov2... # ydotsOjb: the output of findYdotsMM() -# minN: if Ni < minN and have covariates, use the latter instead of # Yi. +# minN: if Ni < minN and have covariates, use the latter instead of Yi. # returns vector of predicted values for testSet predict.ydotsMM = function(ydotsOjb,testSet,minN=0) { haveCovs <- ncol(testSet) > 2 - # use of as.character() is to take advantage of row names, in case of - # future gaps in consecutive user IDs; not implemented yet + # see comment on as.character() above ts1 <- as.character(testSet[,1]) if (!haveCovs) tmp <- ydotsOjb$usrMeans[ts1] else { tmp <- ifelse ( From a7493f4eaf309311922a176b6ffe1d5e559a2280 Mon Sep 17 00:00:00 2001 From: Norm Matloff Date: Sun, 28 May 2017 10:43:03 -0700 Subject: [PATCH 20/30] update man pg, more on filling the gaps --- R/findYdotsMM.R | 5 ++++- man/findYdots.Rd | 5 ++++- 2 files changed, 8 insertions(+), 2 deletions(-) diff --git a/R/findYdotsMM.R b/R/findYdotsMM.R index d275a41..ba25e5a 100644 --- a/R/findYdotsMM.R +++ b/R/findYdotsMM.R @@ -36,7 +36,10 @@ findYdotsMM <- function(ratingsIn,regressYdots=FALSE) { users <- ratingsIn[,1] items <- ratingsIn[,2] # user and item IDs may not be consecutive; even if they are - # consecutive in the original + # consecutive in the original, if we do cross-validation, this + # may not be the case; so switch to character IDs + users <- as.character(users) + items <- as.character(items) ratings <- ratingsIn[,3] nms <- names(ratingsIn) Y.. <- mean(ratings) diff --git a/man/findYdots.Rd b/man/findYdots.Rd index 33430f4..e361bac 100644 --- a/man/findYdots.Rd +++ b/man/findYdots.Rd @@ -16,7 +16,7 @@ capability. Various plots can be displayed. } \usage{ -findYdotsMM(ratingsIn,regressYdots=FALSE,cls=NULL) +findYdotsMM(ratingsIn,regressYdots=FALSE,minN) trainMM(ratingsIn,regressYdots=FALSE,cls=NULL) predict.ydotsMM(ydotsObj,testSet) plot.ydotsMM(ydotsObj,ratingsIn) @@ -31,6 +31,9 @@ predict.ydotsMLE(ydotsObj,testSet) \item{regressYdots}{If TRUE, apply linear regression to the latent factors.} \item{cls}{An R \code{parallel} cluster.} + \item{minN}{If a prediction is to be made on a user with fewer than + this number of ratings and if there are covariates, predict from + the covariates.} \item{ydotsObj}{An object of class \code{'ydotsMM'} or \code{'ydotsMLE'}.} \item{testSet}{A data frame consisting of cases to be predicted. Format is the same as \code{ratingsIn}, except that there is no From f2c85b5bafd41479b811161adbf1a8a866f35f83 Mon Sep 17 00:00:00 2001 From: vchakrab Date: Thu, 1 Jun 2017 20:18:48 -0700 Subject: [PATCH 21/30] hybrid changes --- R/.Rapp.history | 0 R/.Rhistory | 0 R/hybrid.R | 15 +++++++++------ 3 files changed, 9 insertions(+), 6 deletions(-) delete mode 100644 R/.Rapp.history delete mode 100644 R/.Rhistory diff --git a/R/.Rapp.history b/R/.Rapp.history deleted file mode 100644 index e69de29..0000000 diff --git a/R/.Rhistory b/R/.Rhistory deleted file mode 100644 index e69de29..0000000 diff --git a/R/hybrid.R b/R/hybrid.R index 819105a..79a215a 100644 --- a/R/hybrid.R +++ b/R/hybrid.R @@ -8,7 +8,7 @@ # arguments: -# origData: training set, object of class 'usrData' (see file +# datain: training set, object of class 'usrData' (see file # findUsrItmData.R) # newData: data point (just one for now) to be predicted, object of # class 'usrDatum' @@ -33,8 +33,8 @@ hybdpredict.usrData <- function(datain,newData,newItem, } else{ appMatrix = Reco() - train_set = data_memory(unlist(ratingsIn[1]), unlist(ratingIn[2]), - unlist(ratingsIn[3]), index1 = 1) + train_set = data_memory(unlist(datain[1]), unlist(datain[2]), + unlist(datain[3]), index1 = TRUE) appMatrix$train(train_set) #get factors after NMF -- need to use $tune() fix later factors <- appMatrix$output(out_memory(), out_memory()) @@ -57,10 +57,13 @@ hybdpredict.usrData <- function(datain,newData,newItem, #itms <- as.list(as.data.frame(t(oneUsr$itms))) tmp <- match(oneUsr$itms[[1]], newItem) if (all(is.na(tmp))) { - stop('no data on this item') + c(NA, NA) + } + else + { + whichOne <- which(!is.na(tmp)) + c(whichOne,oneUsr$ratings[whichOne]) } - whichOne <- which(!is.na(tmp)) - c(whichOne,oneUsr$ratings[whichOne]) } found <- as.matrix(sapply(basis,checkNewItem)) # found is of dimensions 2 x number of users; From 95e579ec8d6143ae0a8bc29ade3f7fe0a086f5ae Mon Sep 17 00:00:00 2001 From: Pooja-Rajkumar Date: Tue, 6 Jun 2017 14:28:55 -0700 Subject: [PATCH 22/30] added test function --- R/findYdotsMM.R | 69 +++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 69 insertions(+) diff --git a/R/findYdotsMM.R b/R/findYdotsMM.R index ba25e5a..411d5d3 100644 --- a/R/findYdotsMM.R +++ b/R/findYdotsMM.R @@ -113,3 +113,72 @@ checkyd <- function() { print(predict(cout,testset,2)) } +buildMatrix <- function(ratingsIn,NAval=NA){ + # deal with possible factors + dmax <- function(d) { + if (is.factor(d)) return(length(levels(d))) + max(d) + } + users = ratingsIn[,1] + movies = ratingsIn[,2] + rating = ratingsIn[,3] + newMatrix = matrix(NA, + nrow = dmax(users), ncol = dmax(movies)) + for(rows in 1:nrow(ratingsIn)){ + newMatrix[ratingsIn[rows,1],ratingsIn[rows,2]] = ratingsIn[rows,3] + } + return (newMatrix) +} + +predict.nm <- function(approxMatrix, testSet){ + fullMatrix <- buildMatrix(trainSet) # Matrix A (Step 1) + naMatrix <- as.data.frame(which(is.na(fullMatrix) == TRUE, arr.ind = TRUE)) + naMatrix$ratings <- NA + preds <- predict.ydotsMM(approxMatrix,naMatrix) # Step 3 + fullMatrix[which(is.na(fullMatrix))] <- preds + fullMatrix[fullMatrix < 0] <- 0 + fullMatrix[which(is.na(fullMatrix))] <- 0 + res <- nmf(as.matrix(fullMatrix),10) # Step 4 + res +} +# Wrapper for recosystem +trainNM <- function(ratingsIn, trainprop = 0.5,cls = NULL, + rnk = 10, recosystem = FALSE,regressYdots=FALSE){ + require(NMF) + library(NMF) + if(recosystem == TRUE){ + source((paste(system.file(package = 'rectools'), + 'recosys/xvalRecoParallel.R', sep = ""))) + res <- xvalReco(ratingsIn,trainprop,cls,rnk) + } + else { + fullMatrix <- buildMatrix(trainSet) # Matrix A (Step 1) + fullMatrix[which(fullMatrix == 0)] = NA + approxMatrix <- findYdotsMM(trainSet) # Matrix V (Step 2) + } + approxMatrix +} + +# check +check.trainNM <- function(ratingsIn){ + + ## 75% of the sample size + smp_size <- floor(0.25 * nrow(ratingsIn)) + + ## set the seed to make your partition reproductible + set.seed(123) + rownew = nrow(ratingsIn) + trainRow = floor(0.25*rownew) + trainidxs = sample(1:rownew,trainRow) + trainSet = ratingsIn[trainidxs,] + trainSet$trainidxs = trainidxs + testSet = ratingsIn[setdiff(1:nrow(ratingsIn),trainSet$trainidxs),] + trainSet = trainSet[,1:3] + testSet[,3] = NA + trainSet <- trainSet[1:3,] + approxMatrix <- trainNM(trainSet) + res <- predict.nm(approxMatrix,testSet,crossVal = TRUE) + +} + + From 2560cc40f32944dba118b329d33e4bc289d1e19a Mon Sep 17 00:00:00 2001 From: Pooja-Rajkumar Date: Tue, 6 Jun 2017 16:09:46 -0700 Subject: [PATCH 23/30] fix small bug in NAmatrix --- R/findYdotsMM.R | 59 ++++++++++++++++++------------------------------- 1 file changed, 22 insertions(+), 37 deletions(-) diff --git a/R/findYdotsMM.R b/R/findYdotsMM.R index 411d5d3..03a6c91 100644 --- a/R/findYdotsMM.R +++ b/R/findYdotsMM.R @@ -112,6 +112,14 @@ checkyd <- function() { testset$x <- c(5,8) print(predict(cout,testset,2)) } +predict.ydotsMM = function(ydotsObj,testSet) { + testSet$pred = ydotsObj$usrMeans[as.character(testSet[,1])] + + ydotsObj$itmMeans[as.character(testSet[,2])] - ydotsObj$grandMean + if (!is.null(ydotsObj$regObj)) + testSet$pred = testSet$pred + predict(ydotsObj$regObj,testSet[,-(1:2)]) + testSet$pred +} + buildMatrix <- function(ratingsIn,NAval=NA){ # deal with possible factors @@ -122,7 +130,7 @@ buildMatrix <- function(ratingsIn,NAval=NA){ users = ratingsIn[,1] movies = ratingsIn[,2] rating = ratingsIn[,3] - newMatrix = matrix(NA, + newMatrix = matrix(NAval, nrow = dmax(users), ncol = dmax(movies)) for(rows in 1:nrow(ratingsIn)){ newMatrix[ratingsIn[rows,1],ratingsIn[rows,2]] = ratingsIn[rows,3] @@ -130,18 +138,8 @@ buildMatrix <- function(ratingsIn,NAval=NA){ return (newMatrix) } -predict.nm <- function(approxMatrix, testSet){ - fullMatrix <- buildMatrix(trainSet) # Matrix A (Step 1) - naMatrix <- as.data.frame(which(is.na(fullMatrix) == TRUE, arr.ind = TRUE)) - naMatrix$ratings <- NA - preds <- predict.ydotsMM(approxMatrix,naMatrix) # Step 3 - fullMatrix[which(is.na(fullMatrix))] <- preds - fullMatrix[fullMatrix < 0] <- 0 - fullMatrix[which(is.na(fullMatrix))] <- 0 - res <- nmf(as.matrix(fullMatrix),10) # Step 4 - res -} # Wrapper for recosystem + trainNM <- function(ratingsIn, trainprop = 0.5,cls = NULL, rnk = 10, recosystem = FALSE,regressYdots=FALSE){ require(NMF) @@ -152,33 +150,20 @@ trainNM <- function(ratingsIn, trainprop = 0.5,cls = NULL, res <- xvalReco(ratingsIn,trainprop,cls,rnk) } else { - fullMatrix <- buildMatrix(trainSet) # Matrix A (Step 1) - fullMatrix[which(fullMatrix == 0)] = NA - approxMatrix <- findYdotsMM(trainSet) # Matrix V (Step 2) + fullMatrix <- buildMatrix(ratingsIn) # Matrix A (Step 1) + + approxMatrix <- findYdotsMM(ratingsIn) # Matrix V (Step 2) + + naMatrix <- as.data.frame(which(is.na(fullMatrix) == TRUE, arr.ind = TRUE)) + naMatrix$ratings <- NA + preds <- predict.ydotsMM(approxMatrix,naMatrix) # Step 3 + fullMatrix[which(is.na(fullMatrix))] <- preds + fullMatrix[fullMatrix < 0] <- 0 + require(NMF) + res <- nmf(as.matrix(fullMatrix),10) # Step 4 } - approxMatrix + res } -# check -check.trainNM <- function(ratingsIn){ - - ## 75% of the sample size - smp_size <- floor(0.25 * nrow(ratingsIn)) - - ## set the seed to make your partition reproductible - set.seed(123) - rownew = nrow(ratingsIn) - trainRow = floor(0.25*rownew) - trainidxs = sample(1:rownew,trainRow) - trainSet = ratingsIn[trainidxs,] - trainSet$trainidxs = trainidxs - testSet = ratingsIn[setdiff(1:nrow(ratingsIn),trainSet$trainidxs),] - trainSet = trainSet[,1:3] - testSet[,3] = NA - trainSet <- trainSet[1:3,] - approxMatrix <- trainNM(trainSet) - res <- predict.nm(approxMatrix,testSet,crossVal = TRUE) - -} From 2b01a60611d5a5d8869ae069b6379433ff983d7c Mon Sep 17 00:00:00 2001 From: Pooja-Rajkumar Date: Tue, 6 Jun 2017 16:10:25 -0700 Subject: [PATCH 24/30] Update findYdotsMM.R --- R/findYdotsMM.R | 8 +------- 1 file changed, 1 insertion(+), 7 deletions(-) diff --git a/R/findYdotsMM.R b/R/findYdotsMM.R index 03a6c91..7c0d747 100644 --- a/R/findYdotsMM.R +++ b/R/findYdotsMM.R @@ -112,13 +112,7 @@ checkyd <- function() { testset$x <- c(5,8) print(predict(cout,testset,2)) } -predict.ydotsMM = function(ydotsObj,testSet) { - testSet$pred = ydotsObj$usrMeans[as.character(testSet[,1])] + - ydotsObj$itmMeans[as.character(testSet[,2])] - ydotsObj$grandMean - if (!is.null(ydotsObj$regObj)) - testSet$pred = testSet$pred + predict(ydotsObj$regObj,testSet[,-(1:2)]) - testSet$pred -} + buildMatrix <- function(ratingsIn,NAval=NA){ From a46dadacbde057e67112a481db93628ffeb96e38 Mon Sep 17 00:00:00 2001 From: Pooja-Rajkumar Date: Tue, 6 Jun 2017 23:51:25 -0700 Subject: [PATCH 25/30] added predict.nm, getAcc, and getDiff --- R/findYdotsMM.R | 70 +++++++++++++++++++++++++++++++++---------------- 1 file changed, 47 insertions(+), 23 deletions(-) diff --git a/R/findYdotsMM.R b/R/findYdotsMM.R index 7c0d747..de1eec0 100644 --- a/R/findYdotsMM.R +++ b/R/findYdotsMM.R @@ -115,27 +115,46 @@ checkyd <- function() { -buildMatrix <- function(ratingsIn,NAval=NA){ - # deal with possible factors - dmax <- function(d) { - if (is.factor(d)) return(length(levels(d))) - max(d) - } - users = ratingsIn[,1] - movies = ratingsIn[,2] - rating = ratingsIn[,3] - newMatrix = matrix(NAval, - nrow = dmax(users), ncol = dmax(movies)) - for(rows in 1:nrow(ratingsIn)){ - newMatrix[ratingsIn[rows,1],ratingsIn[rows,2]] = ratingsIn[rows,3] - } - return (newMatrix) +predict.NM = function(ydotsObj,testSet) { + testSet$pred = ydotsObj$usrMeans[as.character(testSet[,1])] + + ydotsObj$itmMeans[as.character(testSet[,2])] - ydotsObj$grandMean +# if (!is.null(ydotsObj$regObj)) + # testSet$pred = testSet$pred + predict(ydotsObj$regObj,testSet[,-(1:2)]) + testSet$pred +} + +getDiff <- function(fullMatrix, filledNMF){ + indicies <- which(!is.na(fullMatrix), arr.ind = TRUE) + diff <- matrix(data = NA,nrow = nrow(indicies), 1) + for(i in 1:nrow(indicies)) + for(j in 1:ncol(indicies)){ + diff[i,1] <- round(fullMatrix[indicies[i,1],indicies[j,2]]) - filledNMF[indicies[i,1],indicies[j,2]] + } + + diff } -# Wrapper for recosystem -trainNM <- function(ratingsIn, trainprop = 0.5,cls = NULL, - rnk = 10, recosystem = FALSE,regressYdots=FALSE){ +getAcc <- function(fullMatrix,filledNMF, threshold = 0.5) { + diff <- abs(getDiff(fullMatrix, filledNMF)) + avgDiff <- colMeans(abs(diff), na.rm = TRUE, dims =1) # average difference is 1.19 + # Remove the na's from diff + # Calculate the values less than the threshold + count = 0; + tester <- na.exclude(diff) + for(i in 1:nrow(diff)) + { + if(diff[i] < 0.5 ) + count = count + 1 + } + acc <- count/ nrow(diff) + acc + +} + +train.NM <- function(ratingsIn, trainprop = 0.5,cls = NULL, + rnk = 10, recosystem = FALSE,regressYdots=FALSE) + { require(NMF) library(NMF) if(recosystem == TRUE){ @@ -144,19 +163,24 @@ trainNM <- function(ratingsIn, trainprop = 0.5,cls = NULL, res <- xvalReco(ratingsIn,trainprop,cls,rnk) } else { - fullMatrix <- buildMatrix(ratingsIn) # Matrix A (Step 1) - + origMatrix <- buildMatrix(ratingsIn) # Matrix A (Step 1) + fullMatrix <- origMatrix approxMatrix <- findYdotsMM(ratingsIn) # Matrix V (Step 2) naMatrix <- as.data.frame(which(is.na(fullMatrix) == TRUE, arr.ind = TRUE)) naMatrix$ratings <- NA - preds <- predict.ydotsMM(approxMatrix,naMatrix) # Step 3 + + preds <- predict.NM(approxMatrix,naMatrix) # Step 3 + fullMatrix[which(is.na(fullMatrix))] <- preds fullMatrix[fullMatrix < 0] <- 0 require(NMF) - res <- nmf(as.matrix(fullMatrix),10) # Step 4 + result <- nmf(as.matrix(fullMatrix),10) # Step 4 + filledNMF <- as.matrix(result@fit@W) %*% as.matrix(result@fit@H) + avgd <- getAcc(fullMatrix, filledNMF) # 1.149 for the movielens 100k data + } - res + filledNMF } From 241106838a889477afc97925728490cf28efedad Mon Sep 17 00:00:00 2001 From: Pooja-Rajkumar Date: Tue, 6 Jun 2017 23:54:03 -0700 Subject: [PATCH 26/30] forgot to add buildmatrix --- R/findYdotsMM.R | 17 +++++++++++++++++ 1 file changed, 17 insertions(+) diff --git a/R/findYdotsMM.R b/R/findYdotsMM.R index de1eec0..2db9f98 100644 --- a/R/findYdotsMM.R +++ b/R/findYdotsMM.R @@ -152,6 +152,23 @@ getAcc <- function(fullMatrix,filledNMF, threshold = 0.5) { } +buildMatrix <- function(ratingsIn,NAval=NA){ + # deal with possible factors + dmax <- function(d) { + if (is.factor(d)) return(length(levels(d))) + max(d) + } + users = ratingsIn[,1] + movies = ratingsIn[,2] + rating = ratingsIn[,3] + newMatrix = matrix(NAval, + nrow = dmax(users), ncol = dmax(movies)) + for(rows in 1:nrow(ratingsIn)){ + newMatrix[ratingsIn[rows,1],ratingsIn[rows,2]] = ratingsIn[rows,3] + } + return (newMatrix) +} + train.NM <- function(ratingsIn, trainprop = 0.5,cls = NULL, rnk = 10, recosystem = FALSE,regressYdots=FALSE) { From a45b4d71aad8b7a1c8e72f5c262b1e4372e76779 Mon Sep 17 00:00:00 2001 From: Pooja-Rajkumar Date: Wed, 7 Jun 2017 00:35:53 -0700 Subject: [PATCH 27/30] Update findYdotsMM.R --- R/findYdotsMM.R | 34 +++++++++++++++++----------------- 1 file changed, 17 insertions(+), 17 deletions(-) diff --git a/R/findYdotsMM.R b/R/findYdotsMM.R index 2db9f98..d103dc9 100644 --- a/R/findYdotsMM.R +++ b/R/findYdotsMM.R @@ -152,23 +152,6 @@ getAcc <- function(fullMatrix,filledNMF, threshold = 0.5) { } -buildMatrix <- function(ratingsIn,NAval=NA){ - # deal with possible factors - dmax <- function(d) { - if (is.factor(d)) return(length(levels(d))) - max(d) - } - users = ratingsIn[,1] - movies = ratingsIn[,2] - rating = ratingsIn[,3] - newMatrix = matrix(NAval, - nrow = dmax(users), ncol = dmax(movies)) - for(rows in 1:nrow(ratingsIn)){ - newMatrix[ratingsIn[rows,1],ratingsIn[rows,2]] = ratingsIn[rows,3] - } - return (newMatrix) -} - train.NM <- function(ratingsIn, trainprop = 0.5,cls = NULL, rnk = 10, recosystem = FALSE,regressYdots=FALSE) { @@ -195,10 +178,27 @@ train.NM <- function(ratingsIn, trainprop = 0.5,cls = NULL, result <- nmf(as.matrix(fullMatrix),10) # Step 4 filledNMF <- as.matrix(result@fit@W) %*% as.matrix(result@fit@H) avgd <- getAcc(fullMatrix, filledNMF) # 1.149 for the movielens 100k data + print(avgd) } filledNMF } +buildMatrix <- function(ratingsIn,NAval=NA){ + # deal with possible factors + dmax <- function(d) { + if (is.factor(d)) return(length(levels(d))) + max(d) + } + users = ratingsIn[,1] + movies = ratingsIn[,2] + rating = ratingsIn[,3] + newMatrix = matrix(NAval, + nrow = dmax(users), ncol = dmax(movies)) + for(rows in 1:nrow(ratingsIn)){ + newMatrix[ratingsIn[rows,1],ratingsIn[rows,2]] = ratingsIn[rows,3] + } + return (newMatrix) +} From 5d90ec5fece2d393f635c9b11c5a3c8e8f02e633 Mon Sep 17 00:00:00 2001 From: Norman Matloff Date: Tue, 27 Jun 2017 14:23:46 -0700 Subject: [PATCH 28/30] fixed ivl ex --- man/findYdots.Rd | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/man/findYdots.Rd b/man/findYdots.Rd index e361bac..6f94826 100644 --- a/man/findYdots.Rd +++ b/man/findYdots.Rd @@ -113,6 +113,11 @@ ivl <- ivl[,c(1,2,7,3:6)] library(dummies) dms <- dummy(ivl$dept) dms <- as.data.frame(dms) +# numeric names won't work, so change +dnms <- names(dms) +for (i in 1:length(dnms)) + dnms[i] <- paste('x',dnms[i],sep='') +names(dms) <- dnms dms$dept2 <- NULL ivl$dept <- NULL ivl <- cbind(ivl,dms) From 92c7fa4d4353e61379be55e3319df8d1265f0226 Mon Sep 17 00:00:00 2001 From: Norman Matloff Date: Tue, 27 Jun 2017 17:21:49 -0700 Subject: [PATCH 29/30] fixing MM --- R/findYdotsMM.R | 60 ++++++++++++++++++++++++++++--------------------- 1 file changed, 35 insertions(+), 25 deletions(-) diff --git a/R/findYdotsMM.R b/R/findYdotsMM.R index d103dc9..ca5e50d 100644 --- a/R/findYdotsMM.R +++ b/R/findYdotsMM.R @@ -1,22 +1,18 @@ # modified version of earlier findYdotsMM() etc., May 27, 2017, with -# new approach to use of covariates +# new approach to use of covariates U_ijk for user i, item j; +# NOTE: the covariates will be centered -# in the basic model +# basic model: # -# Y_ij = mu + alpha_i + beta_j + eps +# Y_ij = +# mu + sum_k gamma_k U_ijk + alpha_i + beta_j + eps # -# break down alpha_i: -# -# alpha_i = sum_k phi_i x_ik + gamma_i - # that means -# E(Y_ij | X_ik, k = 1,...,p) = sum_k phi_i x_ik +# E(Y_ij | U_ik, V_jk) = sum_k gamma_k U_ik + sum_k delta_k V_jk -# so the phi can be estimated by lm() without a constant term - -# the same could be done for the items, not yet implemented +# and the coefficients can be estimated via lm() without a const term # arguments: @@ -27,14 +23,15 @@ # S3 class of type "ydots", with components: -# Y..: grand mean (0 if have covariates) -# Yi.: vector of mean ratings for each user -# Y.j: vector of mean ratings for each item -# regOjb: if have covariates, regression output, e.g. coefs +# Y..: grand mean, est of mu +# Yi.: vector of mean ratings for each user, ests. of alpha_i +# Y.j: vector of mean ratings for each item, ests. of betaa_j +# regObj: if have covariates, regression output, e.g. coefs findYdotsMM <- function(ratingsIn,regressYdots=FALSE) { users <- ratingsIn[,1] items <- ratingsIn[,2] + # IMPORTANT NOTE: # user and item IDs may not be consecutive; even if they are # consecutive in the original, if we do cross-validation, this # may not be the case; so switch to character IDs @@ -48,6 +45,10 @@ findYdotsMM <- function(ratingsIn,regressYdots=FALSE) { ydots <- list(grandMean=Y..,usrMeans=Yi.,itmMeans=Y.j) haveCovs <- ncol(ratingsIn) > 3 if (haveCovs) { + # center the covs + tmp <- scale(ratingsIn[,-(1:2)],scale=FALSE) + ratingsIn[,-(1:2)] <- tmp + ydots$covmeans <- attr(tmp,'scaled:center') # regress, no constant term; could do a weighted least squares, # using the Ni, but since the latter is random too, not needed frml <- as.formula(paste(nms[3],'~ .-1')) @@ -74,23 +75,32 @@ trainMM <- findYdotsMM # testSet: data frame in same form as ratingsIn above except that there # is no ratings column; thus covariates, if any, are shifted # leftward one slot, i.e. userID, itemID, cov1, cov2... -# ydotsOjb: the output of findYdotsMM() -# minN: if Ni < minN and have covariates, use the latter instead of Yi. +# ydotsObj: the output of findYdotsMM() +# minN: if Ni < minN and have covariates, use the latter instead of Yi; +# see above # returns vector of predicted values for testSet -predict.ydotsMM = function(ydotsOjb,testSet,minN=0) { +predict.ydotsMM = function(ydotsObj,testSet,minN=0) { haveCovs <- ncol(testSet) > 2 # see comment on as.character() above - ts1 <- as.character(testSet[,1]) - if (!haveCovs) tmp <- ydotsOjb$usrMeans[ts1] else { - tmp <- ifelse ( - ydotsOjb$Ni[ts1] >= minN, - ydotsOjb$usrMeans[ts1], - predict(ydotsOjb$lmout,testSet[,-(1:2),drop=FALSE]) + ts1 <- as.character(testSet[,1]) # user IDs, char form + # tmp will basically consist of the user means, except that in the + # covariate case some will be replaced by predict.lm() values + if (!haveCovs) { + tmp <- ydotsObj$usrMeans[ts1] + } + else { + covs <- testSet[,-(1:2)] + colmeans <- ydotsO + tmp <- + ifelse (ydotsObj$Ni[ts1] >= minN, + ydotsObj$usrMeans[ts1], + predict(ydotsObj$lmout,testSet[,-(1:2),drop=FALSE]) ) } + # so tmp is either estimated mu + alpha or mu + alpha testSet$pred <- tmp + - ydotsOjb$itmMeans[as.character(testSet[,2])] - ydotsOjb$grandMean + ydotsObj$itmMeans[as.character(testSet[,2])] - ydotsObj$grandMean testSet$pred } From 15c046acade2894c48fec43eb0c19691ef610ad3 Mon Sep 17 00:00:00 2001 From: Norm Matloff Date: Tue, 27 Jun 2017 20:53:19 -0700 Subject: [PATCH 30/30] altered model --- R/findYdotsMM.R | 17 +++++++++-------- 1 file changed, 9 insertions(+), 8 deletions(-) diff --git a/R/findYdotsMM.R b/R/findYdotsMM.R index ca5e50d..1bb4aba 100644 --- a/R/findYdotsMM.R +++ b/R/findYdotsMM.R @@ -46,8 +46,8 @@ findYdotsMM <- function(ratingsIn,regressYdots=FALSE) { haveCovs <- ncol(ratingsIn) > 3 if (haveCovs) { # center the covs - tmp <- scale(ratingsIn[,-(1:2)],scale=FALSE) - ratingsIn[,-(1:2)] <- tmp + tmp <- scale(ratingsIn[,-(1:3)],scale=FALSE) + ratingsIn[,-(1:3)] <- tmp ydots$covmeans <- attr(tmp,'scaled:center') # regress, no constant term; could do a weighted least squares, # using the Ni, but since the latter is random too, not needed @@ -85,20 +85,21 @@ predict.ydotsMM = function(ydotsObj,testSet,minN=0) { # see comment on as.character() above ts1 <- as.character(testSet[,1]) # user IDs, char form # tmp will basically consist of the user means, except that in the - # covariate case some will be replaced by predict.lm() values + # covariate case some will be replaced by predict.lm() values + Y.. if (!haveCovs) { tmp <- ydotsObj$usrMeans[ts1] } else { - covs <- testSet[,-(1:2)] - colmeans <- ydotsO + colmeans <- ydotsObj$covmeans + testSet[,-(1:2)] <- + scale(testSet[,-(1:2)],center=colmeans,scale=FALSE) tmp <- ifelse (ydotsObj$Ni[ts1] >= minN, ydotsObj$usrMeans[ts1], - predict(ydotsObj$lmout,testSet[,-(1:2),drop=FALSE]) - ) + predict(ydotsObj$lmout,testSet[,-(1:2),drop=FALSE]) + + ydotsObj$grandMean + ) } - # so tmp is either estimated mu + alpha or mu + alpha testSet$pred <- tmp + ydotsObj$itmMeans[as.character(testSet[,2])] - ydotsObj$grandMean testSet$pred