Skip to contents

This function updates the class means (independent from the other classes).

Usage

update_b(beta, Omega, z, m, xi, Dinv)

Arguments

beta

The matrix of the decision-maker specific coefficient vectors of dimension P_r x N. Set to NA if P_r = 0.

Omega

The matrix of class covariance matrices as columns of dimension P_r*P_r x C. Set to NA if P_r = 0.

z

The vector of the allocation variables of length N. Set to NA if P_r = 0.

m

The vector of class sizes of length C.

xi

The mean vector of length P_r of the normal prior for each b_c. Per default, xi = numeric(P_r).

Dinv

The precision matrix (i.e. the inverse of the covariance matrix) of dimension P_r x P_r of the normal prior for each b_c.

Value

A matrix of updated means for each class in columns.

Details

The following holds independently for each class \(c\). Let \(b_c\) be the mean of class number \(c\). A priori, we assume that \(b_c\) is normally distributed with mean vector \(\xi\) and covariance matrix \(D\). Let \((\beta_n)_{z_n=c}\) be the collection of \(\beta_n\) that are currently allocated to class \(c\), \(m_c\) the class size, and \(\bar{b}_c\) their arithmetic mean. Assuming independence across draws, \((\beta_n)_{z_n=c}\) has a normal likelihood of $$\prod_n \phi(\beta_n \mid b_c,\Omega_c),$$ where the product is over the values \(n\) for which \(z_n=c\) holds. Due to the conjugacy of the prior, the posterior \(\Pr(b_c \mid (\beta_n)_{z_n=c})\) follows a normal distribution with mean $$(D^{-1} + m_c\Omega_c^{-1})^{-1}(D^{-1}\xi + m_c\Omega_c^{-1}\bar{b}_c)$$ and covariance matrix $$(D^{-1} + m_c \Omega_c^{-1})^{-1}.$$

Examples

### N = 100 decider, P_r = 2 random coefficients, and C = 2 latent classes
N <- 100
(b_true <- cbind(c(0,0),c(1,1)))
#>      [,1] [,2]
#> [1,]    0    1
#> [2,]    0    1
Omega <- matrix(c(1,0.3,0.3,0.5,1,-0.3,-0.3,0.8), ncol=2)
z <- c(rep(1,N/2),rep(2,N/2))
m <- as.numeric(table(z))
beta <- sapply(z, function(z) rmvnorm(b_true[,z], matrix(Omega[,z],2,2)))
### prior mean vector and precision matrix (inverse of covariance matrix)
xi <- c(0,0)
Dinv <- diag(2)
### updated class means (in columns)
update_b(beta = beta, Omega = Omega, z = z, m = m, xi = xi, Dinv = Dinv)
#>            [,1]      [,2]
#> [1,]  0.1033711 0.8836673
#> [2,] -0.0073631 0.8993314