0

我正在尝试在加权网络(network_ex)上运行一个有价值的指数随机图模型(ergm)。该网络显示了四个不同组中的个人之间的交互。组之间不能发生交互,因此需要在模型中包含块对角约束。然而,当我包含这个块对角约束时,计算会在蒙特卡洛最大似然估计的特定迭代上卡住至少 24 小时(之后我已经停止了这个过程)。当我删除约束时,计算只需要大约 30 秒。我很困惑,因为我希望块对角约束的计算速度更快,因为可能变化的边数更少。这是一个类似于我的模型的玩具示例:

library(ergm)
library(ergm.count)
library(tergm)

# define response matrix
mat_ex=matrix(
c(0,1,5,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,3,5,12,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,5,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,4,0,0,4,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,7,12,1,5,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,10,3,6,11,6,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,4,0,1,2,6,13,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,5,2,0,11,9,10,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,6,3,0,4,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,6,3,4,8,0,7,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,3,8,1,5,8,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,4,0,0,2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,6,9,11,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,5,3,7,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,5,3,0,9,9,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,10,2,3,0,11,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,5,3,4,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,3,3,0,2,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,4,0,7,0,3,3,1,9,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,5,6,0,7,1,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,5,0,3,0,1,0,1,2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,5,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,3,2,0,0,1,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,1,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,3,1,4,0,0,2,0)
,25,25)

nb=1:25
ID=rep("ID_",25)
names=paste0(ID,nb)
rownames(mat_ex)=names
colnames(mat_ex)=names

# define useful nodes attributes
SeAg=c("F_Adult","F_Adult","M_Adult","M_Adult","M_Sub-Adult","F_Adult","M_Adult","F_Adult","M_Sub-Adult","M_Adult","F_Adult","M_Adult","M_Adult","M_Adult","F_Adult","F_Adult","F_Adult","F_Adult","F_Adult","M_Sub-Adult","M_Adult","M_Adult","M_Adult","F_Adult","M_Adult")
Group=c(1,1,1 ,1,1 ,2, 2, 2 ,2, 2, 2, 2, 3, 3, 3, 3, 3, 4 ,4 ,4, 4, 4, 4,4, 4)
jjj=c(0.9195853,0.5267635,1.0521188,1.4936954,83.8427083,11.3156841,0.6956844,0.2678993,51.5156250,0.7075267,26.5343467,0.7918371,0.3461749,1.5377506,283.7933438,3.2720311,0.9384706,192.6997677,19.0026389,61.6710180,69.0940213,126.9722367,0.5726383,393.5911824,0.7805982)

# define useful edges attributes
mat_iii=matrix(
  c(0.0000000000,0.6413537173,1.7611808824,1.1045365579,0.0610428348,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.6413537173,0.0000000000,0.6214644129,1.0141235877,4.9439358178,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,1.7611808824,0.6214644129,0.0000000000,0.2849059121,0.0013686734,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,1.1045365579,1.0141235877,0.2849059121,0.0000000000,0.0374005963
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0610428348,4.9439358178,0.0013686734,0.0374005963
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.4951618629,0.5858325592,3.8541140580,3.7315606330,0.0009005127
,0.3271637800,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.4951618629,0.0000000000,0.2980769040,0.9923034427,1.4379618200
,3.4747962774,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.5858325592,0.2980769040,0.0000000000,0.2162076270
,0.6839987193,0.0000000000,0.1628661860,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,3.8541140580,0.9923034427,0.2162076270
,0.0000000000,1.2295260710,0.0000000000,1.1504634987,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,3.7315606330,1.4379618200
,0.6839987193,1.2295260710,0.0000000000,1.8073008815,0.5011267262,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0009005127
,3.4747962774,0.0000000000,0.0000000000,1.8073008815,0.0000000000,0.0506384403,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.3271637800,0.0000000000,0.1628661860,1.1504634987,0.5011267262,0.0506384403,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,1.6109005375,4.0166662575,0.4233161826,0.1060254672,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,1.6109005375,0.0000000000,0.3408496603,0.4853374582,0.8949786399,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,4.0166662575,0.3408496603,0.0000000000,0.0698993674,0.6698120875,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.4233161826,0.4853374582,0.0698993674,0.0000000000,1.3822143417
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.1060254672,0.8949786399,0.6698120875,1.3822143417
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.3009462716,2.3082969479,0.1018110048,0.0326341498,0.3828057647
,5.5739028034,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.3009462716,0.0000000000,0.2237564859,0.0323276664,3.8922582745
,0.5733125772,0.6159767242,1.1723116376,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,2.3082969479,0.2237564859,0.0000000000,2.3784384336
,0.7850055818,1.5181554486,0.4208603649,1.7763122536,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.1018110048,0.0323276664,2.3784384336
,0.0000000000,0.1936917389,1.7158463101,0.0000000000,0.0289622499,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0326341498,3.8922582745
,0.7850055818,0.1936917389,0.0000000000,0.8168095460,0.9342629045,0.7772395917,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.3828057647
,0.5733125772,1.5181554486,1.7158463101,0.8168095460,0.0000000000,0.0115138397,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,5.5739028034,0.6159767242,0.4208603649,0.0000000000,0.9342629045,0.0115138397,0.0000000000,1.4325614287
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000,0.0000000000
,0.0000000000,0.0000000000,1.1723116376,1.7763122536,0.0289622499,0.7772395917,0.0000000000,1.4325614287
,0.0000000000)
,25,25)

mat_count_ex=matrix(c(0,1,8,2,2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,4,6,13,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,6,4,0,2,2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,4,0,0,8,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,10,13,2,6,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,16,3,13,13,12,2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,6,0,3,6,9,22,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,6,2,0,14,11,17,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,4,10,3,0,4,3,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,7,3,5,10,0,14,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,4,10,1,6,10,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,9,0,4,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,7,12,16,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,5,4,10,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,8,3,0,13,13,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,14,2,4,0,14,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,5,3,5,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,5,3,1,3,0,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,6,0,7,1,4,4,3,11,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,6,8,0,8,2,0,3,2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,7,2,7,0,2,2,3,9,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,4,5,5,0,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,5,7,0,1,2,0,3,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,2,0,1,2,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,5,2,6,1,0,5,0),25,25)

# define response network
network_ex=as.network(x = mat_ex,
                            directed = TRUE,
                            loops = FALSE,
                            matrix.type = "adjacency",
                            ignore.eval=FALSE,
                            names.eval='weight')
set.vertex.attribute(network_ex, "Sex_Age", as.character(SeAg))
set.vertex.attribute(network_ex, "Group", as.character(Group))
set.vertex.attribute(network_ex, "JJJ", as.numeric(jjj))
set.edge.attribute(network_ex, "III", as.numeric(mat_iii))
set.edge.attribute(network_ex, "Nbobs", as.numeric(mat_count_ex))

# run ergm 
ergm_ex <- ergm(network_ex ~ 
                       sum
                      +nonzero
                      +mutual(form="min",threshold=0)
                      +nodemix("Sex_Age",levels=c(1,2),form="sum")
                      +nodeicov("JJJ",form="sum")
                      +diff("JJJ", pow=1, dir="t-h", sign.action="identity", form ="sum")
                      +edgecov(mat_iii,"III",form="sum") 
                      +edgecov(mat_count_ex,"Nbobs",form="sum")
                      ,response="weight"
                      ,reference=~Poisson
                      ,constraints = ~blockdiag("Group")
                      ,control = control.ergm(MCMC.interval = 1000
                                              ,MCMLE.maxit = 200
                                              ,init.method = 'CD'
                                              ,MCMC.samplesize = 1000
                                              ,MCMC.prop.weights="random"
                                              ,MCMC.burnin=100,seed=123456
                                              ,parallel=4))

这就是我得到的:

最佳有效提案“DistRLE”不能考虑提示“稀疏”。通过 CD-MCMLE 开始对比散度估计:迭代 1,最多 60 次:收敛测试 P 值:0e+00 对数似然提高 0.9066。最多 60 次的第 2 次迭代:收敛检验 P 值:1.5e-248 对数似然提高了 0.4615。最多 60 次的第 3 次迭代:收敛检验 P 值:1.5e-176 对数似然提高了 0.6357。最多 60 次的第 4 次迭代:收敛检验 P 值:6.1e-147 对数似然提高了 0.4542。最多 60 次的第 5 次迭代:收敛检验 P 值:1.1e-92 对数似然提高了 0.2675。最多 60 次的第 6 次迭代:收敛检验 P 值:2.4e-18 对数似然提高了 0.04343。最多 60 次的第 7 次迭代:收敛检验 P 值:1.9e-05 对数似然提高了 0.01517。最多 60 次迭代 8:收敛性检验 P 值:4e-04 对数似然提高了 0.01191。最多 60 次的第 9 次迭代:收敛检验 P 值:6e-03 对数似然提高了 0.009803。最多 60 次的第 10 次迭代:收敛检验 P 值:1.7e-03 对数似然提高了 0.01049。最多 60 次的第 11 次迭代:收敛检验 P 值:6.2e-02 对数似然提高了 0.006888。最多 60 次的第 12 次迭代:收敛检验 P 值:6e-02 对数似然提高了 0.007122。最多 60 次的第 13 次迭代:收敛检验 P 值:2.9e-01 对数似然提高了 0.004809。最多 60 次的第 14 次迭代:收敛检验 P 值:1.5e-02 对数似然提高了 0.00833。最多 60 次的第 15 次迭代:收敛检验 P 值:1.6e-01 对数似然提高了 0.005631。最多 60 次迭代 16:收敛检验 P 值:4。3e-01 对数似然提高了 0.004001。最多 60 次迭代 17:收敛测试 P 值:5.8e-01 检测到收敛。停止。对数似然提高了 0.003348。完成CD。开始 Monte Carlo 最大似然估计 (MCMLE):迭代 1,最多 200 次:优化步长 0.0863421108600662。对数似然提高了 2.006。估计方程不在公差范围内。最多 200 次迭代 2:优化步长 0.0875588561417017。对数似然提高了 2.751。估计方程不在公差范围内。最多 200 次迭代 3:最多 200 次迭代 1:优化步长 0.0863421108600662。对数似然提高了 2.006。估计方程不在公差范围内。最多 200 次迭代 2:优化步长 0.0875588561417017。对数似然提高了 2.751。估计方程不在公差范围内。最多 200 次迭代 3:最多 200 次迭代 1:优化步长 0.0863421108600662。对数似然提高了 2.006。估计方程不在公差范围内。最多 200 次迭代 2:优化步长 0.0875588561417017。对数似然提高了 2.751。估计方程不在公差范围内。最多 200 次迭代 3:

这运行得很快,但不再进步。

4

1 回答 1

0

几点建议:

  • 开发版默认使用自适应 MCMC,并自动执行许多其他操作。我建议在(但不包括)之前删除所有控制参数seed=。(也许seed也是,因为并行处理下的再现性不是很好。)
  • 要更好地了解函数正在做什么,请传递verbose=TRUEorverbose=2或更高。(更高的数字意味着更多的输出,但它们也会减慢速度。)
  • 要查看采样器在做什么,请传递MCMC.runtime.traceplot=TRUE控制参数。
  • 不幸的是,这个采样器效率不高,所以实际上你可能需要运行很长时间。(有一张票:https ://github.com/statnet/ergm.count/issues/5 。)
于 2021-03-31T22:46:37.900 回答