使用ggplot2在R中绘制几个随机试验

Question

我目前正在使用正态分布运行模拟，它模拟事件之间的时间并基于给定数据的分析（与问题无关）。仿真创建这样的：

SimProcess <- function(mu, sigma, T) { 

    ctimes <- c() # Array of arrival times, initially empty 
    t <- rnorm(1,mu, sqrt(sigma)) # Time of next arrival 
    while(t < T) { 
    ctimes <- c(ctimes, t) 
    dt = rnorm(1, mu, sqrt(sigma)) 
    if (dt<0){dt = 0} 
    t <- t + dt # sampling from the dataset 
    } 
    return(ctimes) 
} 

# Create a sample path of one run 
T <- 10 
# arrival times 
arrivals <- SimProcess(mu_t, var_t, T) 

现在我想办几个这样的随机临床试验，然后绘制它们的身影，所以我们可以把它比作给定的数据。这些试验中有10项理想。我试图这样绘图，但不幸的是它不起作用。恐怕我不得不使用reshape2来融合10次试验的数据，因为这些载体的长度都不相同。我用它来试图绘制所有的线条，但它显然不能按照它的方式工作。

x <- c(0, arrivals, T,rep(0,500-length(arrivals))) 
y <- c(0:length(arrivals), length(arrivals),rep(0,500-length(arrivals))) 
plotdataNT = data.frame(x,y) 
p = ggplot(plotdataNT,aes(x,y)) 
plot(x,y,type = 's') 
j = 1 
for (j in 10){ 
    arrivals <- SimProcess(mu_t,var_t,T) 
    x <- c(0, arrivals, T,rep(0,500-length(arrivals))) 
    y <- c(0:length(arrivals), length(arrivals),rep(0,500-length(arrivals))) 
    p = p + geom_step(mapping = aes (x,y)) 
} 

编辑： 最后我figuered出来，因为我用10代替1:10，将无法正常运行，我也提了一些微小的错误。这结束了该解决方案：

arrivals <- SimProcess(mu_t,var_t,T) 
NT <- length(arrivals) 
x <- c(0, arrivals, T,rep(0,correction-length(arrivals))) 
y <- c(0:length(arrivals), length(arrivals),rep(0,correction-length(arrivals))) 
plotdataNT = data.frame(x,y) 
p = ggplot(plotdataNT,aes(x,y)) + geom_step(mapping = aes (x,y)) 
jk = 1 
runs = 25 
colourvec = rainbow(runs) 
for (jk in 1:runs){ 
    arrivals <- SimProcess(mu_t,var_t,T) 
    x <- c(0, arrivals, T,rep(0,correction-length(arrivals))) 
    y <- c(0:length(arrivals), length(arrivals),rep(0,correction-length(arrivals))) 
    newdata = data.frame(x,y) 
    p = p + geom_step(mapping = aes (x,y),newdata,colour = colourvec[jk]) 
} 
p = p + scale_x_continuous(name = "Time in days") + scale_y_continuous(name = "Amount of claims") 
p 

这导致了在几种颜色在一个图中绘制26个的随机样本，它表示具有根据伽马，正常或对数正态分布的随机时间的步骤的过程。下面的答案是我的意思的更清晰的例子。如果有人知道如何以更高效的方式重塑reshape2，我也很乐意知道。

Answer 1

两个解决方案：

for (j in 1:10) { 
    arrivals <- SimProcess(mu_t,var_t,T) 
    x <- c(0, arrivals, T,rep(0,500-length(arrivals))) 
    y <- c(0:length(arrivals), length(arrivals),rep(0,500-length(arrivals))) 
    xy <- data.frame(x,y) 
    p = p + geom_step(data=xy, mapping=aes(x,y)) 
} 
print(p) 

for (j in 1:10) { 
    arrivals <- SimProcess(mu_t,var_t,T) 
    x <- c(0, arrivals, T,rep(0,500-length(arrivals))) 
    y <- c(0:length(arrivals), length(arrivals),rep(0,500-length(arrivals))) 
    xy <- data.frame(x,y) 
    p = p + geom_step(mapping=aes_string(x,y)) 
} 
print(p)