如何将控制权传递给 R 中 for 循环中的上一次迭代
How to pass control to previous iteration in for loop in R
我是 运行 R 中的一个 for 循环(作为模型 I 运行 和 R2jags 的功率分析的一部分)。在某些时候,我想知道我的 MCMC 链是否已经收敛,如果没有,我想跳过循环的迭代。但是,我不想跳到下一次迭代,我希望循环再次从相同的迭代开始。我目前正在使用命令 'next',但这是在跳过迭代。我如何告诉我的 for 循环进行额外的迭代?下面是完整的代码,但我担心的基本上就是这一点:
if(up1 > 1.1 | up2 > 1.1 | up3 > 1.1 | up4 > 1.1)
next
这是完整的代码:
G.vec <- rep(NA, 1000)
#-----model------
model1 <- function(){
# likelihood
for (c in 1:16){
D.diff[c] ~ dnorm(mu, tau)
}
# priors
mu ~ dnorm(0, 10)
tau <- 1/(sd*sd)
sd ~ dunif(0, 10)
}
#---- for loop--------
for (i in 1:1000){
#--- data simulation ---------
# paired data
sim_bf <- rtnorm(16, mean = 0.4949, sd = 0.12, lower = 0.1)
sim_af <- rtnorm(16, mean = 0.3959, sd = 0.12, lower = 0.1)
diff = sim_bf - sim_af
#------setting up jags--------------
data <- list(D.diff = diff)
params <- c("mu", "tau", "sd")
inits <- function(){
list(mu = rnorm(1),
sd = rlnorm(1))
}
#-------run jags----------------
output <- jags(data=data,
# inits=inits,
parameters.to.save=params,
n.iter=1000,
n.burn=100,
n.chains=2,
n.thin=1,
model.file=model1,
progress.bar = "gui")
#-------convergence checker----------------
output.mcmc <- as.mcmc(output)
x <- gelman.diag(output.mcmc)
up1 <- x$psrf[1,2] # Approximate convergence is diagnosed when the upper limit is close to 1
up2 <- x$psrf[2,2]
up3 <- x$psrf[3,2]
up4 <- x$psrf[4,2]
if(up1 > 1.1 | up2 > 1.1 | up3 > 1.1 | up4 > 1.1)
next
# one sided t-test
lo = output$BUGSoutput$summary[2,4]
G.vec[i] <- ifelse((lo < 0), 0, 1)
}
a <- table(G.vec)
G <- a[2]/1000
G
使用 break 语句作为收敛标准。
将下面的行放入 for 循环是没有意义的,因为当收敛满足时它预计 运行。
lo = output$BUGSoutput$summary[2,4]
G.vec[i] <- ifelse((lo < 0), 0, 1)
将它们放在循环之外。
而不是
for(i in 1:1000){...
...
if(...
next
做
i <- 0
while(i <=1000){...
i <- i+1
...
if(...
i <- i-1
break
我是 运行 R 中的一个 for 循环(作为模型 I 运行 和 R2jags 的功率分析的一部分)。在某些时候,我想知道我的 MCMC 链是否已经收敛,如果没有,我想跳过循环的迭代。但是,我不想跳到下一次迭代,我希望循环再次从相同的迭代开始。我目前正在使用命令 'next',但这是在跳过迭代。我如何告诉我的 for 循环进行额外的迭代?下面是完整的代码,但我担心的基本上就是这一点:
if(up1 > 1.1 | up2 > 1.1 | up3 > 1.1 | up4 > 1.1)
next
这是完整的代码:
G.vec <- rep(NA, 1000)
#-----model------
model1 <- function(){
# likelihood
for (c in 1:16){
D.diff[c] ~ dnorm(mu, tau)
}
# priors
mu ~ dnorm(0, 10)
tau <- 1/(sd*sd)
sd ~ dunif(0, 10)
}
#---- for loop--------
for (i in 1:1000){
#--- data simulation ---------
# paired data
sim_bf <- rtnorm(16, mean = 0.4949, sd = 0.12, lower = 0.1)
sim_af <- rtnorm(16, mean = 0.3959, sd = 0.12, lower = 0.1)
diff = sim_bf - sim_af
#------setting up jags--------------
data <- list(D.diff = diff)
params <- c("mu", "tau", "sd")
inits <- function(){
list(mu = rnorm(1),
sd = rlnorm(1))
}
#-------run jags----------------
output <- jags(data=data,
# inits=inits,
parameters.to.save=params,
n.iter=1000,
n.burn=100,
n.chains=2,
n.thin=1,
model.file=model1,
progress.bar = "gui")
#-------convergence checker----------------
output.mcmc <- as.mcmc(output)
x <- gelman.diag(output.mcmc)
up1 <- x$psrf[1,2] # Approximate convergence is diagnosed when the upper limit is close to 1
up2 <- x$psrf[2,2]
up3 <- x$psrf[3,2]
up4 <- x$psrf[4,2]
if(up1 > 1.1 | up2 > 1.1 | up3 > 1.1 | up4 > 1.1)
next
# one sided t-test
lo = output$BUGSoutput$summary[2,4]
G.vec[i] <- ifelse((lo < 0), 0, 1)
}
a <- table(G.vec)
G <- a[2]/1000
G
使用 break 语句作为收敛标准。
将下面的行放入 for 循环是没有意义的,因为当收敛满足时它预计 运行。
lo = output$BUGSoutput$summary[2,4] G.vec[i] <- ifelse((lo < 0), 0, 1)
将它们放在循环之外。
而不是
for(i in 1:1000){...
...
if(...
next
做
i <- 0
while(i <=1000){...
i <- i+1
...
if(...
i <- i-1
break