The R package mfair implements the methods based on the paper MFAI: A
scalable Bayesian matrix factorization approach to leveraging auxiliary
information. MFAI integrates
gradient boosted trees in the probabilistic matrix factorization
framework to leverage auxiliary information effectively and adaptively.
For a quick start, you can install the development version of mfair
from GitHub with:
# install.packages("devtools")
devtools::install_github("YangLabHKUST/mfair")For more illustration and examples, you can alternatively use:
# install.packages("devtools")
devtools::install_github("YangLabHKUST/mfair", build_vignettes = TRUE)to build vignettes simultaneously. Please note that it can take a few more minutes.
- This is a basic example which shows you how to solve a common problem:
set.seed(20230306)
library(mfair)
# Simulate data
# Set the data dimension and rank
N <- 100
M <- 100
K_true <- 2L
# Set the proportion of variance explained (PVE)
PVE_Z <- 0.8
PVE_Y <- 0.5
# Generate auxiliary information X
X1 <- runif(N, min = -10, max = 10)
X2 <- runif(N, min = -10, max = 10)
X <- cbind(X1, X2)
# F(X)
FX1 <- X1 / 2 - X2
FX2 <- (X1^2 - X2^2 + 2 * X1 * X2) / 10
FX <- cbind(FX1, FX2)
# Generate the factor matrix Z (= F(X) + noise)
sig1_sq <- var(FX1) * (1 / PVE_Z - 1)
Z1 <- FX1 + rnorm(n = N, mean = 0, sd = sqrt(sig1_sq))
sig2_sq <- var(FX2) * (1 / PVE_Z - 1)
Z2 <- FX2 + rnorm(n = N, mean = 0, sd = sqrt(sig2_sq))
Z <- cbind(Z1, Z2)
# Generate the loading matrix W
W <- matrix(rnorm(M * K_true), nrow = M, ncol = K_true)
# Generate the main data matrix Y_obs (= Y + noise)
Y <- Z %*% t(W)
Y_var <- var(as.vector(Y))
epsilon_sq <- Y_var * (1 / PVE_Y - 1)
Y_obs <- Y + matrix(
rnorm(N * M,
mean = 0,
sd = sqrt(epsilon_sq)
),
nrow = N, ncol = M
)
# Create MFAIR object
mfairObject <- createMFAIR(Y_obs, as.data.frame(X), K_max = K_true)
#> The main data matrix Y is completely observed!
#> The main data matrix Y has been centered with mean = 0.147726471347656!
# Fit the MFAI model
mfairObject <- fitGreedy(mfairObject, sf_para = list(verbose_loop = FALSE))
#> Set K_max = 2!
#> Initialize the parameters of Factor 1......
#> After 2 iterations Stage 1 ends!
#> After 77 iterations Stage 2 ends!
#> Factor 1 retained!
#> Initialize the parameters of Factor 2......
#> After 2 iterations Stage 1 ends!
#> After 76 iterations Stage 2 ends!
#> Factor 2 retained!
# Prediction based on the low-rank approximation
Y_hat <- predict(mfairObject)
#> The main data matrix Y has no missing entries!
# Root-mean-square-error
sqrt(mean((Y_obs - Y_hat)^2))
#> [1] 12.23344
# Predicted/true matrix variance ratio
var(as.vector(Y_hat)) / var(as.vector(Y_obs))
#> [1] 0.4714952
# Prediction/noise variance ratio
var(as.vector(Y_hat)) / var(as.vector(Y_obs - Y_hat))
#> [1] 0.9871629mfaircan also handle the matrix with missing entries:
# Split the data into the training set and test set
n_all <- N * M
training_ratio <- 0.5
train_set <- sample(1:n_all, n_all * training_ratio, replace = FALSE)
Y_train <- Y_test <- Y_obs
Y_train[-train_set] <- NA
Y_test[train_set] <- NA
# Create MFAIR object
mfairObject <- createMFAIR(Y_train, as.data.frame(X), K_max = K_true)
#> The main data matrix Y is partially observed!
#> The main data matrix Y has been centered with mean = 0.187847085351627!
# Fit the MFAI model
mfairObject <- fitGreedy(mfairObject, sf_para = list(verbose_loop = FALSE))
#> Set K_max = 2!
#> Initialize the parameters of Factor 1......
#> After 2 iterations Stage 1 ends!
#> After 97 iterations Stage 2 ends!
#> Factor 1 retained!
#> Initialize the parameters of Factor 2......
#> After 2 iterations Stage 1 ends!
#> After 82 iterations Stage 2 ends!
#> Factor 2 retained!
# Prediction based on the low-rank approximation
Y_hat <- predict(mfairObject)
# Root-mean-square-error
sqrt(mean((Y_test - Y_hat)^2, na.rm = TRUE))
#> [1] 13.08502
# Predicted/true matrix variance ratio
var(as.vector(Y_hat), na.rm = TRUE) / var(as.vector(Y_obs), na.rm = TRUE)
#> [1] 0.4078598
# Prediction/noise variance ratio
var(as.vector(Y_hat), na.rm = TRUE) / var(as.vector(Y_obs - Y_hat), na.rm = TRUE)
#> [1] 0.7989475- Empirically, the backfitting algorithm can further improve the performance:
# Refine the MFAI model with the backfitting algorithm
mfairObject <- fitBack(mfairObject,
verbose_bf_inner = FALSE,
sf_para = list(verbose_sf = FALSE, verbose_loop = FALSE)
)
#> Iteration: 1, relative difference of model parameters: 0.2678141.
#> Iteration: 2, relative difference of model parameters: 0.03957596.
#> Iteration: 3, relative difference of model parameters: 0.08902799.
#> Iteration: 4, relative difference of model parameters: 0.02089378.
#> Iteration: 5, relative difference of model parameters: 0.001688755.
# Prediction based on the low-rank approximation
Y_hat <- predict(mfairObject)
# Root-mean-square-error
sqrt(mean((Y_test - Y_hat)^2, na.rm = TRUE))
#> [1] 13.03505
# Predicted/true matrix variance ratio
var(as.vector(Y_hat), na.rm = TRUE) / var(as.vector(Y_obs), na.rm = TRUE)
#> [1] 0.4259078
# Prediction/noise variance ratio
var(as.vector(Y_hat), na.rm = TRUE) / var(as.vector(Y_obs - Y_hat), na.rm = TRUE)
#> [1] 0.8400624vignette("ml100k")- Explore the vignette illustrating the spatial and temporal dynamics of gene regulation among brain tissues:
vignette("neocortex")- For more documentation and examples, please visit our package website.
If you find the mfair package or any of the source code in this
repository useful for your work, please cite:
Wang, Z., Zhang, F., Zheng, C., Hu, X., Cai, M., and Yang, C. (2023). MFAI: A scalable Bayesian matrix factorization approach to leveraging auxiliary information. arXiv preprint arXiv:2303.02566. URL: https://doi.org/10.48550/arXiv.2303.02566.
The package is developed by Zhiwei Wang ([email protected]).
Please feel free to contact Zhiwei Wang ([email protected]), Prof. Mingxuan Cai ([email protected]), or Prof. Can Yang ([email protected]) with any inquiries.