R - 在循环内循环答案

【问题标题】：R - Looping within a loopR - 在循环内循环
【发布时间】：2019-02-28 18:49:10
【问题描述】：

我使用下面的代码来计算一个函数的对数似然，其中 x 参数 (x1, x2, x3, x4, x5) 是在 0 和 1 之间均匀分布的随机数，但该函数还依赖于t 的值（它是一个向量）。

library(expm)
t <- seq(10,1000,by=5)

x1 <- runif(1,0,1)
x2 <- runif(1,0,1)
x3 <- runif(1,0,1)
x4 <- runif(1,0,1)
x5 <- runif(1,0,1)

p <- matrix(c(1,0,0), nrow=1, ncol=3, byrow=TRUE)

q <- matrix(c(x1,x2,x3), nrow=3, ncol=1, byrow=TRUE)

likely <- vector()

for (i in 1:length(t)) {
likely[i] <- p%*%expm(matrix(c(-(x4+x1)*t[i],x4*t[i],0,0,-(x5+x2)*t[i],x5*t[i],0,0,-x3*t[i]), nrow=3, ncol=3, byrow=TRUE))%*%q
}

log.likely <- vector()

for (i in 1:length(t)) {
log.likely[i] <- log(likely[i])
}

L3 <- -(sum(log.likely))
L3

我的问题是：对于 x 参数的不同值，我该如何运行此代码 100 次？我已经在代码中有一个循环来遍历 t 的不同值。但现在我需要遍历 x1、x2、x3、x4 和 x5 的不同值。任何帮助将不胜感激！

【问题讨论】：

创建一个包含您要检查的所有参数组合的数据框或矩阵。 expand.grid 对此很有用，例如 params = expand.grid(x1 = c(0.1, 0.5, 0.72), x2 = seq(0, 1, by = 0.2), ...)。然后创建一个遍历1:nrow(params) 的循环，运行你已经为每一行拥有的代码，并存储结果。
如果您将现有代码转换为以x1-x5 作为参数的函数，这可能是最简单的。
@Gregor 感谢您的回复！我现在有这样的东西：params <- expand.grid(x1 = runif(10,0,1), x2 = runif(10,0,1), x3 = runif(10,0,1), x4 = runif(10,0,1), x5 = runif(10,0,1)) 但我不确定如何在循环中存储矩阵。例如，这不起作用：for (k in 1:nrow(params)) { q[k] <- matrix(c(x1,x2,x3), nrow=3, ncol=1, byrow=TRUE) }
q是长度列表=nrow(params)吗？这可能行得通。

标签： r loops

【解决方案1】：

我不确定这是否正是您想要的，但这将允许您多次运行您的代码，并将结果存储在包含所有 x 值的列表中以及这些值的对数可能性。

log.likelyhood <- function(){
  require(expm)

  t <- seq(10,1000,by=5)

  x1 <- runif(1,0,1)
  x2 <- runif(1,0,1)
  x3 <- runif(1,0,1)
  x4 <- runif(1,0,1)
  x5 <- runif(1,0,1)

  p <- matrix(c(1,0,0), nrow=1, ncol=3, byrow=TRUE)

  q <- matrix(c(x1,x2,x3), nrow=3, ncol=1, byrow=TRUE)

  likely <- vector()

  for (i in 1:length(t)) {
    likely[i] <- p%*%expm(matrix(c(-(x4+x1)*t[i],x4*t[i],0,0,-(x5+x2)*t[i],x5*t[i],0,0,-x3*t[i]), nrow=3, ncol=3, byrow=TRUE))%*%q
  }

  log.likely <- vector()

  for (i in 1:length(t)) {
    log.likely[i] <- log(likely[i])
  }

  L3 <- -(sum(log.likely))

  return(list(xvals = c(x1, x2, x3, x4, x5), L3 = L3))
}

results <- replicate(100, log.likelyhood(), simplify = FALSE)

【讨论】：

这正是我要找的，谢谢！它在我提供的时间向量上运行完美。但是，我的实时向量有大量元素，所以我仍在等待生成结果。我不确定这需要多长时间，但再次感谢您的回答！