hierMeanReg | R Documentation |
fits a hierarchical normal model of the form E[y_{ij}] = μ_{j} + β_{1}x_{i1}+…+β_{p}x_{ip}
hierMeanReg( design, priorTau, priorPsi, priorVar, priorBeta = NULL, steps = 1000, startValue = NULL, randomSeed = NULL )
design |
a list with elements y = response vector, group = grouping vector, x = matrix of covariates or NULL if there are no covariates |
priorTau |
a list with elements tau0 and v0 |
priorPsi |
a list with elements psi0 and eta0 |
priorVar |
a list with elements s0 and kappa0 |
priorBeta |
a list with elements b0 and bMat or NULL if x is NULL |
steps |
the number of Gibbs sampling steps to take |
startValue |
a list with possible elements tau, psi, mu, sigmasq and beta. tau, psi and sigmasq must all be scalars. mu and beta must be vectors with as many elements as there are groups and covariates respectively |
randomSeed |
a random seed for the random number generator |
A data frame with variables:
tau |
Samples from the posterior distribution of tau |
psi |
Samples from the posterior distribution of psi |
mu |
Samples from the posterior distribution of mu |
beta |
Samples from the posterior distribution of beta if there are any covariates |
sigmaSq |
Samples from the posterior distribution of σ^2 |
sigma |
Samples from the posterior distribution of sigma |
priorTau = list(tau0 = 0, v0 = 1000) priorPsi = list(psi0 = 500, eta0 = 1) priorVar = list(s0 = 500, kappa0 = 1) priorBeta = list(b0 = c(0,0), bMat = matrix(c(1000,100,100,1000), ncol = 2)) data(hiermeanRegTest.df) data.df = hiermeanRegTest.df design = list(y = data.df$y, group = data.df$group, x = as.matrix(data.df[,3:4])) r=hierMeanReg(design, priorTau, priorPsi, priorVar, priorBeta) oldPar = par(mfrow = c(3,3)) plot(density(r$tau)) plot(density(r$psi)) plot(density(r$mu.1)) plot(density(r$mu.2)) plot(density(r$mu.3)) plot(density(r$beta.1)) plot(density(r$beta.2)) plot(density(r$sigmaSq)) par(oldPar) ## example with no covariates priorTau = list(tau0 = 0, v0 = 1000) priorPsi = list(psi0 = 500, eta0 = 1) priorVar = list(s0 = 500, kappa0 = 1) data(hiermeanRegTest.df) data.df = hiermeanRegTest.df design = list(y = data.df$y, group = data.df$group, x = NULL) r=hierMeanReg(design, priorTau, priorPsi, priorVar) oldPar = par(mfrow = c(3,2)) plot(density(r$tau)) plot(density(r$psi)) plot(density(r$mu.1)) plot(density(r$mu.2)) plot(density(r$mu.3)) plot(density(r$sigmaSq)) par(oldPar)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.