R获取概率分布

Question

我有一个关系式:y=a+b+c

我有a，b和c的平均值和标准差，我想通过蒙特卡罗模拟得到y的概率分布。

有没有什么函数、包或者简单的方法可以让我这样做呢？

Answer 1

我假设你的输入a，b和c是正态分布的，因为你说你可以用均值和标准差来定义它们。如果是这种情况，您可以在没有任何特殊包的情况下快速完成此操作。

 mu.a=33
 mu.b=32
 mu.c=13
 sigma.a=22
 sigma.b=22
 sigma.c=222

n= a.large.number=10^5
a=rnorm(n,mu.a,sigma.a)
b=rnorm(n,mu.b,sigma.b)
c=rnorm(n,mu.c,sigma.c)
y=a+b+c
plot(density(y))
mean(y)
sd(y)

确保了解我们对y、a、b和c所做的所有假设。如果你想做一些更复杂的事情，比如计算出y的均值的抽样方差，然后多次执行这个过程，收集均值并绘制出来。

mysimfun=function(n,mu,sigma,stat.you.want='mean') 
   #  mu is length 3 and sigma is too.

{
n= a.large.number=10^5
    a=rnorm(n,mu[1],sigma[1])
    b=rnorm(n,mu[2],sigma[2])
    c=rnorm(n,mu[3],sigma[3])
    y=a+b+c
    plot(density(y))


return(ifelse(stat.you.want=='mean',mean(y),sd(y))
}


mu=c(mu.a,my.b,mu.c)
sigma=c(sigma.a,sigma.b,sigma.c)
mi=rep(NA,100)

然后在某种循环中运行它。

for(i in 1:100) {mi[i]=mysimfun(10,mu,sigma,stat.you.want='mean') }

par(mfrow=c(2,1)
hist(mi)
plot(density(mi))

mean(mi)
sd(mi)

Answer 2

有两种方法: bootstrapping，我认为这就是你所说的MonteCarlo，或者如果你对这个理论比从经验分布中构建估计更感兴趣，'distr‘包及其朋友’distrSim‘和"distrTEst“。

require(boot)
ax <- rnorm(100); bx<-runif(100); cx<- rexp(100)
dat <- data.frame(ax=ax,bx=bx,cx=cx)

boot(dat, function(d){ with(d, mean(ax+bx+cx) )}, R=1000,  sim="parametric")
boot(dat, function(d){ with(d, sd(ax+bx+cx) )}, R=1000,  sim="parametric")