Description Usage Arguments Value Examples
unimod
is mainly modeling the measures for combining summaries from multiple sources. It returns the estimates and posterior distributions of parameters of interest fitted with data.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 |
y |
A n-by-1 vector, consisting of the measures of n observations. |
s |
A n-by-1 vector, consisting of the standard error of the measures for n observations. |
with_cov |
if TRUE, include covariates in the model. |
X |
A n-by-p matrix of p covariates including intercept, which is all 1 for the first column. |
seed |
The seed for random number generation. The default is generated from 1 to the maximum integer supported by R on the machine. Even if multiple chains are used, only one seed is needed, with other chains having seeds derived from that of the first chain to avoid dependent samples.The default is 9999. |
iter |
a number to specify the iteration times of Stan. The default is 3000. |
verbose |
TRUE or FALSE: flag indicating whether to print intermediate output from Stan on the console, which might be helpful for model debugging. |
warmup |
a number of warmup iterations to be excluded when computing the summaries. The default is 1000. |
thin |
A positive integer specifying the period for saving samples. The default value is 3. |
adapt_delta |
A parameters to control the sampler’s behavior. The default is 0.8. |
max_treedepth |
A parameters to control the sampler’s behavior in For algorithm NUTS. |
chains |
A positive integer specifying the number of Markov chains. The default is 4. |
a list of information of the fitted univariate Bayestian hierarchical model including
pop_par
the estimates and credible intervals for the population-level parameters
ind_par
the estimates and credible intervals for the individual-level parameters
poster_dist
posterior distributions of all parameters fitted in the model
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 | library(rstan)
options(mc.cores = parallel::detectCores())
rstan_options(auto_write = TRUE)
library(MASS)
set.seed(99999)
### define data generation functions-with predictors
m <- 50
X=t(rbind(rep(1,m),rnorm(m,0,1),rbinom(m,1,0.2)))
beta=t(rbind(c(4,3,5),c(1,1,0)))
data.reg <- function(tautheta=3,tausigma=1,
r1=0.3,r2=0.3, sigmas=1,...){
# r1 - individual level
# r2 - population level
n <- dim(X)[1]
theta = vector()
sigma2 = vector()
y = vector()
s2 = vector()
Sigma = list()
temp = matrix(rep(NA),ncol = 2, nrow = n)
indiv = matrix(rep(NA),ncol = 2, nrow = n)
# generate individual true mean and variance from population
mu <- X%*%beta
varmat <- matrix(c(tautheta^2,r2*tausigma*tautheta,r2*tausigma*tautheta,tausigma^2), nrow=2)
for (k in 1:n) temp[k,] <- MASS::mvrnorm(1,mu[k,],varmat)
theta <- temp[,1]
sigma2 <- (exp(temp[,2]))^2
# generate observed value from true individual
for (i in 1:n){
Sigma[[i]] <- c(sigma2[i],r1*sqrt(sigma2[i])*sigmas,r1*sqrt(sigma2[i])*sigmas,sigmas^2)
indiv[i,] <- mvrnorm(1, temp[i,],matrix(Sigma[[i]], nrow = 2))
}
y <- indiv[,1]
s <- exp(indiv[,2])
dataset <- data.frame(cbind(y, s, X))
return(list(data = dataset,
para = c(tautheta=tautheta, tausigma=tausigma, r1=r1, r2=r2, sigmas=sigmas)))
}
data <- data.reg(r1=0,r2=0)$data
res <- unimod(data$y, data$s, with_cov = TRUE, X=X,
seed = 9999, iter = 3000, chain = 3, verbose = FALSE,
warmup= 1000, thin=3)$pop_par
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.