#-------------
# y ~ Gaussian
#-------------
# model assumption:
theta <- c(1, 2, 3, 4, 1.5) # coefficients and sigma2 = 1.5
#----------------
# Data Simulation
#----------------
n <- 40
X <- data.frame(x1=rnorm(n), # continuous variable
x2=sample(1:3, n, replace=TRUE)) # categorical variable
Xx2_1 <- ifelse(X$x2 == '2', 1, 0)
Xx2_2 <- ifelse(X$x2 == '3', 1, 0)
X$x2 <- as.factor(X$x2)
eta <- theta[1] + theta[2] * X$x1 + theta[3] * Xx2_1 + theta[4] * Xx2_2
mu <- gaussian()$linkinv(eta)
y <- rnorm(n, mu, sd = sqrt(theta[5]))
#----------------
# MAP estimations
#----------------
Lambda <- inv.prior.cov(X, lambda = c(0.1, 0.5), family = "gaussian")
fit <- MAP.estimation(y, X, family = "gaussian", Lambda)
class(fit)
#-------------------------
# Summary of MAP estimates
#-------------------------
summary(fit)
sumfit <- summary(fit, cur_mat = TRUE)
sumfit$estimate
sumfit$logLikPost
sumfit$dispersion
sumfit$CI
class(sumfit)
Run the code above in your browser using DataLab