在 r 中模拟复合泊松过程

Question

我正在尝试在 r 中模拟复合泊松过程。该过程由 $ \sum_{j=1}^{N_t} Y_j $ 定义，其中 $Y_n$ 是 i.i.d 序列独立 $N(0,1) $ values 和 $N_t$ 是参数为 $1$ 的泊松过程。我试图在 r 中模拟这个但没有运气。我有一个算法来计算这个如下：模拟从 0 到 T 的 cPp：

发起：$k=0$

重复 $\sum_{i=1}^k T_i < T$

设 $k = k+1$

模拟 $T_k \sim exp(\lambda)$（在我的例子中 $\lambda = 1$）

Simulate $Y_k \sim N(0,1)$（这只是一个特例，我希望能够将其更改为任何分布）

轨迹由 $X_t = \sum_{j=1}^{N_t} Y_j $ 给出，其中 $N(t) = sup(k : \ sum_{i=1}^k T_i \leq t )$

谁能帮我在 r 中模拟一下，这样我就可以绘制出这个过程？我已经尝试过，但无法完成。

Answer 1

使用 cumsum 作为确定时间 N_t 和 X_t 的累积总和。此说明性代码指定要模拟的次数 n，模拟 n.t 中的时间和 x 中的值，并（显示它所做的）绘制轨迹。

n <- 1e2
n.t <- cumsum(rexp(n))
x <- c(0,cumsum(rnorm(n)))
plot(stepfun(n.t, x), xlab="t", ylab="X")

这个算法，因为它依赖于 low-level 优化函数，速度很快：我测试它的 six-year-old 系统每秒将生成超过三百万（时间，值）对。

这对于模拟来说通常已经足够好了，但它并不能完全满足问题，该问题要求生成一个到时间 T 的模拟。我们可以利用前面的代码，但解决方案有点棘手。它计算在时间 T 之前泊松过程中将发生多少次的合理上限。它生成 inter-arrival 次。这被包裹在一个循环中，该循环将在（罕见的）事件中重复该过程，但实际未达到时间 T。

额外的复杂度不会改变渐近计算时间。

T <- 1e2            # Specify the end time
T.max <- 0          # Last time encountered
n.t <- numeric(0)   # Inter-arrival times
while (T.max < T) {
  #
  # Estimate how many random values to generate before exceeding T.
  #
  T.remaining <- T - T.max
  n <- ceiling(T.remaining + 3*sqrt(T.remaining))
  #
  # Continue the Poisson process.
  #
  n.new <- rexp(n)
  n.t <- c(n.t, n.new)
  T.max <- T.max + sum(n.new)
}
#
# Sum the inter-arrival times and cut them off after time T.
#
n.t <- cumsum(n.t)
n.t <- n.t[n.t <= T]
#
# Generate the iid random values and accumulate their sums.
#
x <- c(0,cumsum(rnorm(length(n.t))))
#
# Display the result.
#
plot(stepfun(n.t, x), xlab="t", ylab="X", sub=paste("n =", length(n.t)))

在 r 中模拟复合泊松过程

Simulate Compound poisson process in r

simulation

r