M
那么,如何用拟合出来的分布做Simulation呢,就要用到一个知识点,那就是可以先生成服从0到1之间均匀分布的随机数,再Map到拟合分段函数的CDF $F(x)$上,然后计算该概率下的quantile就是最后的随机数。
分段函数的CDF $F(x)$是由每段函数原本的CDF通过改变截距和形状来的,那么我们就可以把它还原到原来函数的CDF,用R的功能计算出quantile来
#生成1000个符合0-1之间均匀分布的随机数
init_numebrs = runif(1000)
simulation = function(init){
#这两个break就是CDF的分段点
break1= 0.2
break2 = 0.2+0.5
#这里就是把改变形状截距的CDF还原回去,然后计算quantile
if(init<break1){
output = qgamma(init/modelCollections[[1]]['coefficient'][[1]],shape=modelCollections[[1]]\$estimate["shape"], rate=modelCollections[[1]]\$estimate["rate"])
}else if(init>break2){
output = qgamma((init-break2)/modelCollections[[3]]['coefficient'][[1]]+pgamma(120,modelCollections[[3]]\$estimate["shape"],modelCollections[[3]]\$estimate["rate"]),shape=modelCollections[[3]]\$estimate["shape"], rate=modelCollections[[3]]\$estimate["rate"])
}else{
output = qweibull((init-break1)/modelCollections[[2]]['coefficient'][[1]]+pweibull(50,modelCollections[[2]]\$estimate["shape"],modelCollections[[2]]\$estimate["scale"]),shape=modelCollections[[2]]\$estimate["shape"], scale=modelCollections[[2]]\$estimate["scale"])
}
}
simulation_result=lapply(init_numebrs,simulation)
dt = data.frame("data"=unlist(simulation_result))
p = ggplot(data = data.frame(dt), aes(x = data)) + geom_histogram(binwidth=10,fill="#69b3a2", color="#e9ecef", alpha=0.3)
拟合出来的分布大概是以下形状