Skip to content

Instantly share code, notes, and snippets.

@RottenFruits

RottenFruits/Bayesian_Probabilistic_Matrix_Factorization.ipynb

Created August 16, 2018 11:27
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save RottenFruits/f0aeb144ea50dfd689362187ef725966 to your computer and use it in GitHub Desktop.
Save RottenFruits/f0aeb144ea50dfd689362187ef725966 to your computer and use it in GitHub Desktop.
Download ZIP
Display the source blob
Display the rendered blob
Raw
Bayesian_Probabilistic_Matrix_Factorization.ipynb
{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# module"
]
},
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"Main.BayesianProbabilisticMatrixFactorization"
]
},
"execution_count": 1,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"module BayesianProbabilisticMatrixFactorization\n",
"\n",
"using LinearAlgebra\n",
"using Distributions\n",
"\n",
"mutable struct BayesianProbabilisticMatrixFactorizationModel\n",
" D::Int64\n",
" alpha::Float64\n",
" beta_0::Float64\n",
" mu_0::Float64\n",
" nu_0::Float64\n",
" W_0::Array\n",
" N::Int64\n",
" M::Int64\n",
" T::Int64\n",
" R::Array\n",
" U::Array\n",
" V::Array\n",
"end\n",
"\n",
"#fit関数\n",
"function fit(model::BayesianProbabilisticMatrixFactorizationModel)\n",
" model.U = zeros((model.T, model.N, model.D)) #ユーザー潜在変数初期化\n",
" model.V = zeros((model.T, model.M, model.D)) #アイテム潜在変数初期化\n",
" \n",
" model.U[1, :, :] = rand(model.N, model.D) #初期値セット\n",
" model.V[1, :, :] = rand(model.M, model.D) #初期値セット\n",
" \n",
" T = model.T\n",
" N = model.N\n",
" M = model.M\n",
" \n",
" for t in 1:(T-1) #1には初期値セットしているので1引く\n",
" #u matrix\n",
" parameter_u = hyper_prameter_posterior(model, model.U[t, :, :])\n",
" mu_u, Lambda_u = hyper_parameter_sampling(parameter_u[1], parameter_u[2], parameter_u[3], parameter_u[4])\n",
"\n",
" #v matrix\n",
" parameter_v = hyper_prameter_posterior(model, model.V[t, :, :])\n",
" mu_v, Lambda_v = hyper_parameter_sampling(parameter_v[1], parameter_v[2], parameter_v[3], parameter_v[4]) \n",
"\n",
" for n in 1:N\n",
" mu_i_asterisk, Lambda_i_asterisk_inv = u_matrix_posterior(model, model.V[t, :, :], mu_u, Lambda_u, n)\n",
" model.U[t + 1, n, :] = latent_matrix_sampling(mu_i_asterisk, Lambda_i_asterisk_inv)\n",
" end\n",
"\n",
" for m in 1:M\n",
" mu_i_asterisk, Lambda_i_asterisk_inv = v_matrix_posterior(model, model.U[t + 1, :, :], mu_v, Lambda_v, m)\n",
" model.V[t + 1, m, :] = latent_matrix_sampling(mu_i_asterisk, Lambda_i_asterisk_inv)\n",
" end\n",
"\n",
" end\n",
"end\n",
"\n",
"#ハイパーパラメータの事後分布\n",
"function hyper_prameter_posterior(model::BayesianProbabilisticMatrixFactorizationModel, latent_matrix)\n",
" N = model.N\n",
" beta_0 = model.beta_0\n",
" mu_0 = model.mu_0\n",
" nu_0 = model.nu_0\n",
" W_0 = model.W_0\n",
" \n",
" latent_matrix_bar = (1 / N) * sum(latent_matrix[:, :], dims =1) \n",
" S_bar = (1 / N) .* sum([latent_matrix[i, :]' .* latent_matrix[i, :] for i = 1:N], dims = 1)[1] \n",
" \n",
" mu_0_asterisk = ((beta_0 * mu_0 .+ N * latent_matrix_bar) / (beta_0 + N))[1, :]\n",
" beta_0_asterisk = beta_0 + N\n",
" nu_0_asterisk = nu_0 + N\n",
" \n",
" W_0_asterisk = inv(inv(W_0) .+ N .* S_bar .+ ((beta_0 * N) / (beta_0 + N)) .* (mu_0 .- latent_matrix_bar)' * (mu_0 .- latent_matrix_bar))\n",
" W_0_asterisk = Array(Symmetric(W_0_asterisk))\n",
" \n",
" return mu_0_asterisk, beta_0_asterisk, nu_0_asterisk, W_0_asterisk \n",
"end\n",
"\n",
"#ハイパーパラメータのサンプリング\n",
"function hyper_parameter_sampling(mu_0_asterisk, beta_0_asterisk, nu_0_asterisk, W_0_asterisk)\n",
" Lambda = rand(Wishart(nu_0_asterisk, W_0_asterisk), 1)[1]\n",
" mu = rand(MvNormal(mu_0_asterisk, Array(Symmetric(inv(beta_0_asterisk * Lambda)))), 1)\n",
" \n",
" return mu, Lambda\n",
"end\n",
"\n",
"#ユーザ潜在行列の事後分布\n",
"function u_matrix_posterior(model::BayesianProbabilisticMatrixFactorizationModel, latent_matrix, mu_u, Lambda_u, n)\n",
" R = model.R\n",
" alpha = model.alpha\n",
" \n",
" r = R[R[:, 1] .== n - 1, :] #ユーザーの評価取り出す\n",
" nonzero_ix = (R[R[:, 1] .== n - 1, :][:, 2] .+ 1) #ユーザーの評価つけているアイテムのインデクス\n",
" \n",
" if nonzero_ix == [] # 評価がない場合は適当な値を返す\n",
" mu_i_asterisk = zeros(size(latent_matrix)[2])\n",
" tmp_matrix = rand(size(latent_matrix)[2], size(latent_matrix)[2])\n",
" Lambda_i_asterisk = Float64.(Int.(inv(tmp_matrix) * tmp_matrix .>= 0.9)) #逆行列かけて型変換をして単位行列に変換\n",
" \n",
" return mu_i_asterisk, Lambda_i_asterisk\n",
" end\n",
"\n",
" Lambda_i_asterisk = Lambda_u .+ alpha * sum([latent_matrix[i, :]' .* latent_matrix[i, :] for i in nonzero_ix], dims = 1)[1]\n",
" Lambda_i_asterisk_inv = Array(Symmetric(inv(Lambda_i_asterisk)))\n",
" \n",
" mu_i_asterisk = (inv(Lambda_i_asterisk) * ((alpha * sum([latent_matrix[i, :] .* r[r[:, 2] .== i - 1, :][:, 3] for i in nonzero_ix], dims = 1)[1]') + (Lambda_u * mu_u)')')[:, 1]\n",
" \n",
" return mu_i_asterisk, Lambda_i_asterisk_inv\n",
"end\n",
"\n",
"#アイテム潜在行列の事後分布\n",
"function v_matrix_posterior(model::BayesianProbabilisticMatrixFactorizationModel, latent_matrix, mu_v, Lambda_v, m)\n",
" R = model.R\n",
" alpha = model.alpha\n",
" \n",
" r = R[R[:, 2] .== m - 1, :] #評価しているユーザーを取り出す\n",
" nonzero_ix = (R[R[:, 2] .== m - 1, :][:, 1] .+ 1) #アイテムに評価つけているユーザーのインデクス\n",
" \n",
" if nonzero_ix == [] # 評価がない場合は適当な値を返す\n",
" mu_i_asterisk = zeros(size(latent_matrix)[2])\n",
" tmp_matrix = rand(size(latent_matrix)[2], size(latent_matrix)[2])\n",
" Lambda_i_asterisk = Float64.(Int.(inv(tmp_matrix) * tmp_matrix .>= 0.9)) #逆行列かけて型変換をして単位行列に変換\n",
" \n",
" return mu_i_asterisk, Lambda_i_asterisk\n",
" end\n",
" \n",
" Lambda_i_asterisk = Lambda_v .+ alpha * sum([latent_matrix[i, :]' .* latent_matrix[i, :] for i in nonzero_ix], dims = 1)[1]\n",
" Lambda_i_asterisk_inv = Array(Symmetric(inv(Lambda_i_asterisk)))\n",
"\n",
" mu_i_asterisk = (inv(Lambda_i_asterisk) * ((alpha * sum([latent_matrix[i, :] .* r[r[:, 1] .== i - 1, :][:, 3] for i in nonzero_ix], dims = 1)[1]') + (Lambda_v * mu_v)')')[:, 1]\n",
" \n",
" return mu_i_asterisk, Lambda_i_asterisk_inv\n",
"end\n",
"\n",
"#潜在行列のサンプリング\n",
"function latent_matrix_sampling(mu_i_asterisk, Lambda_i_asterisk_inv)\n",
" return rand(MvNormal(mu_i_asterisk, Lambda_i_asterisk_inv), 1)\n",
"end\n",
"\n",
"#全予測\n",
"function predict_all(model::BayesianProbabilisticMatrixFactorizationModel, burn_in)\n",
" T = model.T\n",
" \n",
" sum([model.U[i, :, :] * model.V[i,:, :]' for i = burn_in:T], dims = 1)[1] / (T - burn_in )\n",
"end\n",
"\n",
"#全予測データを使って新規データ予測\n",
"function predict(new, R_p)\n",
" r = zeros(size(new)[1])\n",
" \n",
" for i in 1:size(new)[1]\n",
" r[i] = R_p[new[i, 1] + 1, new[i, 2] + 1]\n",
" end\n",
" \n",
" return r\n",
"end\n",
"\n",
"end"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# experiment"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"using DataFrames\n",
"using CSV\n",
"using Random\n",
"using StatsModels\n",
"using Statistics"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"20001×3 Array{Int64,2}:\n",
" 483 37 4\n",
" 6 594 2\n",
" 805 482 4\n",
" 649 1038 3\n",
" 402 865 4\n",
" 906 276 5\n",
" 853 602 4\n",
" 245 474 4\n",
" 871 120 4\n",
" 560 173 4\n",
" 329 20 5\n",
" 14 147 3\n",
" 124 708 3\n",
" ⋮ \n",
" 540 28 2\n",
" 330 681 5\n",
" 289 379 3\n",
" 478 1243 3\n",
" 659 82 3\n",
" 348 712 3\n",
" 488 268 3\n",
" 42 472 3\n",
" 541 55 5\n",
" 587 70 4\n",
" 607 1182 1\n",
" 560 659 3"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"#データ読み込み\n",
"df = CSV.read(\"ml-100k/u.data\", header = false, delim = '\\t')\n",
"\n",
"#配列化\n",
"df = convert(Array{Int64}, df[:, 1:3])\n",
"\n",
"#オフセット idの最小値を0にする\n",
"df[:, 1, :] = df[:, 1, :] .- 1\n",
"df[:, 2, :] = df[:, 2, :] .- 1\n",
"\n",
"#ユニークユーザー、ユニークアイテム\n",
"user = length(unique(df[:, 1]))\n",
"item = length(unique(df[:, 2]))\n",
"\n",
"#シャッフル\n",
"df = df[shuffle(1:end), :]\n",
"\n",
"#学習データとテストデータ分割\n",
"N = size(df)[1]\n",
"train_size = Int64(N * 0.8)\n",
"train_df = df[1:train_size, :]\n",
"test_df = df[train_size:N, :]"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"209.223477 seconds (212.96 M allocations: 304.652 GiB, 9.71% gc time)\n"
]
},
{
"data": {
"text/plain": [
"943×1682 Array{Float64,2}:\n",
" 3.88961 2.97194 3.07 4.22053 … 3.81322 2.56164 4.05672 3.43431\n",
" 4.10209 3.36757 2.58199 3.69592 3.88245 2.72459 4.22869 3.5284 \n",
" 2.96916 2.54889 3.40775 3.22184 3.32122 2.06453 3.60701 3.06467\n",
" 5.09577 4.35659 3.65003 4.22654 4.93761 3.37646 5.56712 4.52131\n",
" 3.64042 3.10722 2.79575 3.66518 3.56167 2.39591 3.80533 3.08185\n",
" 3.44175 2.61876 1.57762 3.5659 … 3.37799 2.20563 3.52767 3.08162\n",
" 4.31396 3.74416 3.58239 4.0239 4.34943 2.96584 4.84158 4.06383\n",
" 4.19003 3.26543 3.15086 3.90335 3.89826 2.77172 4.20783 3.46531\n",
" 4.70969 4.23153 3.95679 3.70391 4.50848 3.13257 5.11352 4.12165\n",
" 4.38406 3.80762 3.37166 3.8115 4.27812 2.91659 4.68036 3.88306\n",
" 3.81955 3.11471 2.98878 3.48148 … 3.92098 2.60121 4.23381 3.505 \n",
" 4.47758 3.60567 3.65307 4.06238 4.23289 2.79175 4.66678 3.73438\n",
" 3.41243 2.78145 2.61979 3.89528 3.42182 2.42914 3.76735 3.21126\n",
" ⋮ ⋱ ⋮ \n",
" 4.03603 3.21275 2.46239 3.74218 3.92312 2.71765 4.16286 3.5473 \n",
" 2.84408 2.32009 1.84603 3.19849 2.89436 1.8586 2.95987 2.59075\n",
" 3.87825 3.2606 3.52144 3.63157 3.9019 2.67303 4.33769 3.54072\n",
" 4.15088 3.92932 3.03493 2.88425 3.73289 2.80742 4.36023 3.40462\n",
" 4.34645 3.54277 3.56106 3.80493 … 4.34818 2.85803 4.74809 3.92728\n",
" 3.91352 3.37489 1.57292 3.3139 3.40164 2.59572 3.74117 3.09782\n",
" 3.7279 3.51558 3.57921 3.44146 3.78975 2.80221 4.3786 3.63753\n",
" 5.11007 4.34774 4.10742 4.42836 4.90928 3.30353 5.3619 4.33199\n",
" 4.12695 3.67028 2.54014 3.03583 3.42825 2.49317 3.88157 3.16798\n",
" 4.72479 3.86559 3.12431 4.12423 … 4.47072 2.92389 4.85833 3.91478\n",
" 4.60918 3.93191 3.08954 3.8279 4.25564 2.84945 4.67992 3.78768\n",
" 3.80317 3.48886 3.37123 3.81054 3.51688 2.32469 3.78156 3.1147 "
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"#時間計測\n",
"function main()\n",
" #hyper hyper parameter\n",
" beta_0 = 2\n",
" mu_0 = 0\n",
" D = 10 #分解後の潜在変数の次元\n",
" nu_0 = D \n",
" tmp_matrix = rand(D, D)\n",
" W_0 = Float64.(Int.(inv(tmp_matrix) * tmp_matrix .>= 0.9)) #逆行列かけて型変換をして単位行列に変換\n",
" alpha = 2\n",
" T = 50 #サンプリング回数\n",
"\n",
" N = user #ユーザー数\n",
" M = item #アイテム数\n",
"\n",
" R = df#今回は簡易的な実験のため学習と精度検証用のデータはわけない\n",
"\n",
" #学習\n",
" BPMF = BayesianProbabilisticMatrixFactorization.BayesianProbabilisticMatrixFactorizationModel(D, alpha, beta_0, mu_0, nu_0,\n",
" W_0, N, M, T, R, [], [])\n",
" BayesianProbabilisticMatrixFactorization.fit(BPMF)\n",
"\n",
" #予測\n",
" R_P = BayesianProbabilisticMatrixFactorization.predict_all(BPMF, 10)\n",
" \n",
" return R_P\n",
"end\n",
"\n",
"@time R_P = main()"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.7190330856378937"
]
},
"execution_count": 7,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"sqrt(mean((df[:, 3] - BayesianProbabilisticMatrixFactorization.predict(df, R_P)) .^ 2))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# 精度検証用データを分けた場合の精度"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"0.9337503598743961"
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"#hyper hyper parameter\n",
"beta_0 = 2\n",
"mu_0 = 0\n",
"D = 10 #分解後の潜在変数の次元\n",
"nu_0 = D \n",
"tmp_matrix = rand(D, D)\n",
"W_0 = Float64.(Int.(inv(tmp_matrix) * tmp_matrix .>= 0.9)) #逆行列かけて型変換をして単位行列に変換\n",
"alpha = 2\n",
"T = 50 #サンプリング回数\n",
"\n",
"N = user #ユーザー数\n",
"M = item #アイテム数\n",
"\n",
"R = train_df\n",
"\n",
"\n",
"#学習\n",
"BPMF = BayesianProbabilisticMatrixFactorization.BayesianProbabilisticMatrixFactorizationModel(D, alpha, beta_0, mu_0, nu_0,\n",
" W_0, N, M, T, R, [], [])\n",
"BayesianProbabilisticMatrixFactorization.fit(BPMF)\n",
"\n",
"#予測\n",
"res = BayesianProbabilisticMatrixFactorization.predict_all(BPMF, 10)\n",
"\n",
"#rmse\n",
"sqrt(mean((test_df[:, 3] - BayesianProbabilisticMatrixFactorization.predict(test_df, res)) .^ 2))"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Julia 1.0.0",
"language": "julia",
"name": "julia-1.0"
},
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "1.0.0"
}
},
"nbformat": 4,
"nbformat_minor": 1
}
Display the source blob
Display the rendered blob
Raw
bpmf_python.ipynb
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment