在r中模拟复合泊松过程

Question

我在模拟r中的复合泊松过程。该过程由$ \ sum_ {j = 1}^{N_t} Y_j $定义，其中$ Y_n $是i.i.d序列无关的$ N（0,1）$值，$ N_t $是参数为$ 1 $的泊松过程。我试图在没有运气的情况下模拟这个。我有一个算法来计算这个如下： Simutale的CPP从0到T：

启动：$ K = 0 $

重复而$ \ sum_ {I = 1} ^ķT_i < T $

设置$ k = k + 1 $

模拟$ T_k \ SIM EXP（\拉姆达）$（在我的情况$ \拉姆达= $ 1）

模拟$ Y_K \ SIM N（0 ，1）$（这只是一个特殊情况，我希望能够将其更改为任何分配）

该轨迹由下式给出：其中$ 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") 

这种算法，因为它依赖于低级别的优化功能，速度快：我测试了它在六十岁系统将生成超过三百万（时间，值）对每秒。

---

这通常用于模拟不够好，但它并不完全符合的问题，这要求产生一个模拟出时间T.我们可以利用上面的代码，但解决的办法是有点麻烦。它计算出在时间T之前泊松过程中会发生多少次的合理上限。它会产生到达时间间隔。这被封装在一个循环中，该循环将在实际上没有达到时间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)))