非常非常感谢！我尝试一下

再请问一下，模型做完之后，怎样画出随机效应的 BLUE ?  以及怎样得出μi～ (μ, σ2μ), βi～ (β, σ2β), εij～ (0, σ2)这三个正态分布的具体的参数来？拜托了

什么意思，我只是按照模型simulate了一个数据。如果还要估计参数的话：

非常感谢，解决了困扰我很久的问题。假定有 n 个个体被观测, 第 i 个个体被观测 mi次, 第 j 次观测结果为
yij= μi+ xijβi+ εij， i = 1, 2, · · · , n; j = 1, 2, · · · , mi. 其中 μi～ (μ, σ2μ), βi～ (β, σ2β), εij～ (0, σ2)
依然以这个为例，比如我按照这个程序做出一个随机系数模型来,固定效应：y=458.67*x-814.51，随机效应： μi～ (μ, 54456), βi～ (β, 15336), εij～ (0, 188)，假定他们都是正态分布，而且我的结果里面 i = 1, 2, · · · , 18; j = 1, 2, · · · , 21. 那么我想请问的是：您上面simulate了一个数据，我想按照我做出来的这个模型随便拿出一个x值然后来simulate y的值，simulate 1000次以上，这样plot（y）应该就可以看出预测值y的分布范围，但是μ和β的均值也没给出，我也不知道这该怎么设定，这样的话程序该怎么写，所以还请您多多指导，不胜感激！