I know usually, we do not use Bayesian Information criterion(BIC) for model selection if we have a linear mixed model (problems involve like the sample size in the linear mixed model is not well defined etc.). But I was trying to see a "rough" BIC result about the linear mixed model, so I was running some numerical experiments about it.

I kind of want to repeat the situations in paper, and compare my result with the result in this paper but with a different method.

Suppose there are $N$ subjects under study, with subject $i$ contribution $n_i$ observations, for $i =1,...,N.$. And let $y_{ij}$ denote a response variable for subject $i$ at observation $j$. Let $x_{ij}$ denote a $p\times 1$ vector of predictors, and let $z_{ij}$ denote a $q\times 1$ vector of predictors. In general, the linear mixed effects model is $$y_i=X_i\alpha+Z_i\beta_i+\epsilon_i$$ where $y_i=(y_{i1},...,y_{in_i})^T$, $X_i=(x_{i1}^T,...,x_{in_i}^T)^T$, $Z_i=(z_{i1}^T,...,z_{in_i}^T)^T$, $\alpha$ is a $p \times 1$ vector of unknown population parameters, $\beta_i$ is a $q \times 1$ vector of unknown subject-specific random effects with $\beta_i \sim N(0, D)$ and the ellements of the residual vector, $\epsilon_i$, are $N(0, \sigma^2I)$.

In particular, we consider the case where $q=4, N=200$, with 8 observations for each subject. The covariates $x_{ij}=(x_{ij1},...,x_{ij4})^T$ are simulated by fixing $x_{ij1}=1$ and then generating $x_{ijk} \sim \text{Uniform}(-2,2)$ for $k=2,3,4.$ We let $z_{ij}=(z_{ij1},...,z_{ij4})^T=x_{ij}$ and choose $\alpha=(1,1,1,1)^T$ and $\sigma^2=1$ in the model with $\beta_i=(\beta_{i1},...,\beta_{i4})^T \sim N(0, D)$, where $D=\begin{pmatrix} 9 & 4.8 & 0.6 & 0 \\ 4.8 & 4 & 1 & 0 \\ 0.6 & 1 & 1 & 0 \\ 0 & 0 & 0 & 0 \end{pmatrix}$.

To simplify the question, we only compare the model $$y_{ij}=\alpha_1+\alpha_2x_{ij2}+\alpha_3x_{ij3}+\alpha_4x_{ij4}+\beta_{i1}+\beta_{i2}z_{ij2}+\beta_{i3}z_{ij3}+\beta_{i4}z_{ij4}+\epsilon_{ij}.$$ with $$y_{ij}=\alpha_1+\alpha_2x_{ij2}+\alpha_3x_{ij3}+\alpha_4x_{ij4}+\beta_{i1}+\beta_{i2}z_{ij2}+\beta_{i3}z_{ij3}+\epsilon_{ij}.$$ What I was doing is generate $y_i \sim N(X_i\alpha, Z_iDZ_i^T+\sigma^2I)$, then generate $X_i$ and $Z_i$, then do the linear regression to see the BIC.

But I am a little confused about what I did. For example, if I just regress $y$ on those $X$'s and $Z$'s, is that the same compare with a non-mixed model? i.e. how can we see the random effect is actually random here?

Also, by forcing $D_{44}=0$, we actually force $\beta_4=0$, i.e. the second model shuld be better than the first one, but the result is actually highly depend on the data, I do not have a lower BIC for the second model always.

Any ideas will be very helpful!


  • $\begingroup$ I cannot understand why we need BIC, given $ \mathrm {BIC} =\ln(n)k-2\ln(\hat L)$ $\endgroup$ – user158565 Dec 2 '18 at 20:15
  • $\begingroup$ @user158565 Just want to see how BIC performs in the mixed model. $\endgroup$ – Nan Dec 2 '18 at 20:25
  • $\begingroup$ In your situation, what is the deference between studying BIC and studying $\ln(\hat L)$ given $n, k_1, k_2=k_1-4$ are constant? $\endgroup$ – user158565 Dec 2 '18 at 20:28
  • $\begingroup$ @user158565 Sorry I did not quite get your point. I just using R to get the log likelihood first then get the BIC. $\endgroup$ – Nan Dec 2 '18 at 20:30
  • $\begingroup$ So just study log-likelihood. To see if LRT can reject the null hypotheses that the fourth random effect does not exist at the given alpha level. $\endgroup$ – user158565 Dec 2 '18 at 20:51

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.