get_moments_cpp_eco: Calculate variational moments during the updates (only for...

In zhenkewu/lotR: Latent class analysis with Observed Trees in R (lotR)

Calculate variational moments during the updates (only for node u)

Description

update only selected moments that need update when iterating over u (except for rmat)

Usage

get_moments_cpp_eco(
  leaves_u,
  E_beta,
  E_beta_sq,
  E_eta,
  E_eta_sq,
  prob,
  prob_gamma,
  mu_gamma,
  sigma_gamma,
  mu_alpha,
  Sigma_alpha,
  anc,
  cardanc
)

Arguments

leaves_u	the leaf descendant node ids for node u
E_beta, E_beta_sq, E_eta, E_eta_sq	moment updates produced by get_moments_cpp
prob	variational probabilities for s_u; length p
prob_gamma	should be fixed: c(1,rep(0,p-1))
mu_gamma	variational Gaussian means (for s_u=1 component) for J*K logit(class-specific response probabilities); (J,K,p) array; In R, we used a list of p (J,K) matrices
sigma_gamma	variational Gaussian variances (for s_u=1 component) for J*K logit(class-specific response probabilities); (J,K,p) array
mu_alpha	variational Gaussian mean vectors (for s_u=1 component) - this is a p by K-1 matrix; in R, we used a list of p vectors (each of length K-1)
Sigma_alpha	variational Gaussian variances (for s_u=1 component) this is an array of dimension (K-1, K-1, p); in R, we used a list of p matrices, each of dimension K-1 by K-1.
anc	a list of pL vectors, each vector has the node ids of the ancestors; lengths may differ. The ancestors include the node concerned.
cardanc	a numeric vector of length pL; integers. The number of ancestors for each leaf node

Details

(one-node version of get_moments_cpp)

Value

a List

return List::create(Named("E_beta")=E_beta, Named("E_beta_sq")=E_beta_sq, Named("E_eta")=E_eta, Named("E_eta_sq")=E_eta_sq);

Examples

# illustrates the calculation of the various moments:

prob = 0.2
mu1 = 2
mu0 = 0

sigma1 = 3
sigma0 = 1

N = 10000
s = rbinom(N,1,prob)
y <- rep(NA,length(s))
for (i in seq_along(s)){
      y[i] <- rnorm(1,mu1,sqrt(sigma1))*s[i] +
        rnorm(1,mu0,sqrt(sigma0))*(1-s[i])
}

var(s*y)
prob*(sigma1+(1-prob)*mu1^2)

mean(y^2)
prob*(sigma1+mu1^2) + (1-prob)*sigma0

mean(s*y^2)
prob*(sigma1+mu1^2)

mean(s*y)
prob*mu1

zhenkewu/lotR documentation built on April 24, 2022, 2:36 a.m.