fosr: MCMC Sampling Algorithm for the Function-on-Scalars...

In drkowal/fosr: Function-on-Scalars Regression for High Dimensional Data

Description

Runs the MCMC for the function-on-scalars regression model based on an FDLM-type expansion. Here we assume the factor regression has independent errors, which allows for subject-specific random effects, as well as some additional default conditions.

Usage

fosr(
  Y,
  tau,
  X = NULL,
  K = NULL,
  nsave = 1000,
  nburn = 1000,
  nskip = 3,
  mcmc_params = list("beta", "fk", "alpha", "sigma_e", "sigma_g", "trsigma"),
  computeDIC = TRUE
)

Arguments

Y	the n x m data observation matrix, where n is the number of subjects and m is the number of observation points (NAs allowed)
tau	the m x d matrix of coordinates of observation points
X	the n x p matrix of predictors; if NULL, only include an intercept
K	the number of factors; if NULL, use SVD-based proportion of variability explained
nsave	number of MCMC iterations to record
nburn	number of MCMC iterations to discard (burin-in)
nskip	number of MCMC iterations to skip between saving iterations, i.e., save every (nskip + 1)th draw
mcmc_params	named list of parameters for which we store the MCMC output; must be one or more of "beta" (factors) "fk" (loading curves) "alpha" (regression coefficients) "sigma_e" (observation error SD) "sigma_g" (random effects SD) "Yhat" (fitted values) "trsigma" (the sum of trace Sigma_i) "sigma_delta_k" (random effects gamma specific SD)
computeDIC	logical; if TRUE, compute the deviance information criterion DIC and the effective number of parameters p_d

Value

A named list of the nsave MCMC samples for the parameters named in mcmc_params

Note

If nm is large, then storing all posterior samples for Yhat, which is nsave x n x M, may be inefficient

Examples

# Simulate some data:
sim_data = simulate_fosr(n = 100, m = 20, p_0 = 100, p_1 = 5)

# Data:
Y = sim_data$Y; X = sim_data$X; tau = sim_data$tau

# Dimensions:
n = nrow(Y); m = ncol(Y); p = ncol(X)

# Run the FOSR:
out = fosr(Y = Y, tau = tau, X = X, K = 6, mcmc_params = list("fk", "alpha", "Yhat"))

# Plot a posterior summary of a regression function, say j = 3:
j = 3; post_alpha_tilde_j = get_post_alpha_tilde(out$fk, out$alpha[,j,])
plot_curve(post_alpha_tilde_j, tau = tau)
# Add the true curve:
lines(tau, sim_data$alpha_tilde_true[,j], lwd=6, col='green', lty=6)

# Plot the loading curves:
plot_flc(out$fk, tau = tau)

# Plot the fitted values for a random subject:
i = sample(1:n, 1)
plot_fitted(y = Y[i,], mu = colMeans(out$Yhat[,i,]),
            postY = out$Yhat[,i,], y_true = sim_data$Y_true[i,], t01 = tau)

drkowal/fosr documentation built on Oct. 2, 2020, 11:20 a.m.