拓端tecdat|R语言连续时间马尔可夫链模拟案例 Markov Chains

原文链接：http://tecdat.cn/?p=4182

 

案例

一个加油站有一个加油桩，没有空间供车辆等待（如果车辆到达，加油桩被占用，它就会离开）。车辆到达加油站的速率服从泊松过程λ=3/20每分钟，其中75％是汽车，25％是摩托车。加油时间可以用一个指数随机变量建模，平均汽车8分钟，摩托车3分钟，服务速率为汽车μC= 1 / 8和摩托车μ= 1 / 3 每分钟。

因此，我们可以通过将这些概率乘以每个状态下的车辆数量来计算系统中的平均车辆数量。

1.  
   # 到达率
2.  
   lambda <- 3/20
3.  
   # 服务速率（汽车，摩托车）
4.  
   mu <- c(1/8, 1/3)
5.  
   # 汽车的概率
6.  
   p <- 0.75
7.  
    
8.  
   #理论解析
9.  
   A <- matrix(c(1, mu[1], 0,
10.  
    1, -lambda, (1-p)*lambda,
11.  
     
12.  
    N_average_theor
13.  
    #> [1] 0.5031056

现在，我们将模拟系统并验证

1.  
   optio<-
2.  
    
3.  
   seize("pump", amount=1) %>%
4.  
   timeout(function() rexp(1, mu[1])) %>%
5.  
   release("pump", amount=1)
6.  
    

为了区分汽车和摩托车，我们可以在获取资源后定义一个分支来选择合适的服务时间。

这option.3相当于option.1性能。 例如，

1.  
   opti2 <- function(t) {
2.  
    
3.  
   seize("pump", amount=1) %>%
4.  
   branch(function() sample(c(1, 2), 1, prob=c(p, 1-p)), c(T, T),
5.  
   trajectory("car")
6.  
   timeout(function() rexp(1, mu[2]))) %>

但是此选项增加了不必要的运算，因为需要额外调用R函数来选择分支，因此会降低性能。更好的选择是直接在timeout函数内部选择服务时间。

1.  
   optio3 <- function(t) {
2.  
   vehicle <- trajectory() %>%
3.  
   seize("pump", amount=1) %>%
4.  
   if (runif(1) < p) rexp(1, mu[1]) # 汽车
5.  
   else rexp(1, mu[2]) # 摩托车
6.  
   }) %>%

这option.3等效option.1于性能。但是，我们得出了相同的结果。例如，

1.  
    
2.  
   # 使用率+理论值
3.  
   plot(get_mon_resources(gas.station), "usage", "pump", items="system") +
4.  
   geom_hline(yintercept=N_average_theor)

 

这些是一些表现的结果：

1.  
    
2.  
   t <- 1000/lambda
3.  
   tm <- microbenchmark(option.1(t),
4.  
    
5.  
    
6.  
    
7.  
   autoplot(tm) +
8.  
   scale_y_log10(breaks=function(limits) pretty(limits, 5)) +

---

最受欢迎的见解

1.用R语言模拟混合制排队随机服务排队系统

2.R语言中使用排队论预测等待时间

3.R语言中实现马尔可夫链蒙特卡罗MCMC模型

4.R语言中的马尔科夫机制转换(Markov regime switching)模型

5.matlab贝叶斯隐马尔可夫hmm模型

6.用R语言模拟混合制排队随机服务排队系统

7.Python基于粒子群优化的投资组合优化

8.R语言马尔可夫转换模型研究交通伤亡人数事故预测

9.用机器学习识别不断变化的股市状况——隐马尔可夫模型的应用