如何用r语言采用蒙特卡罗模拟5000次计算出EWMA控制图的平均运行链长ARL

附上我参考的文献及我自己写的代码，根据参考文献写的，但是运行结果反过来了 ，求大神援助帮我改对代码就可以，事关毕业论文，谢谢啦，重赏！！！1
#指数分布，参数为R
#正态分布0.9975对应 的分位点为2.81,即ARL0取值为200时对应的L值为2.81
#lambda=1,样本数为200，每个样本容量为7，即可控下的数据位200，第201个数据开始发生偏移。
lambda=0.25
L=2.81
shift=c(1,1.1,1.2,1.3,1.4,1.5,1.6,1.8,2,2.5,3,3.5,4)#标准差偏移度i
l=length(shift) #13个不同的偏移度
n=100#j
newsigma=matrix(,n,l)
rl=matrix(1300,n,l)
for(i in 1:l){
       for(j in 1:n){
       m=matrix(rexp(1400,1),200,7)  #生成200个样本容量为7，参数lambda为1的指数分布
       #参数R
       R=c()
       Rsigma=c()
       Rz=c()
       sigma=sd(m)
       for(k in 1:nrow(m)){
       R[k]=max(m[k,])-min(m[k,])#求每个样本极差
       Rz[k]=R[k]/sigma           #求每个样本参数Rz
       }
       Rbar=mean(R)
       Rt2=mean(Rz)
       Rt3=sd(Rz)
       LCL=Rbar-(L*Rt3*Rbar*(lambda/(2-lambda))^0.5)/Rt2  #控制下限
       UCL=Rbar+(L*Rt3*Rbar*(lambda/(2-lambda))^0.5)/Rt2  #控制上限
       sdw=(Rt3*Rbar*(lambda/(2-lambda))^0.5)/Rt2   #计算过程方差
     newsigma[j,i]=shift*sdw      #发生偏移后的标准差
     ewma=c()
     num=1
     ewma[1]=Rbar
     for(k in 2:nrow(m)){
    ewma[k]=lambda*newsigma[j,i]+(1-lambda)*ewma[k-1] #发生偏移后的统计量
     if (ewma[k]>LCL&&ewma[k]<UCL){                      #判断是否在区间范围内
        num=num+1           #计算当出现第一个不在控制区间时已经落在控制上下限的总数
        }
        else{break}

  }
  rl[j,i]=num #求每一次的运行链长
}

}
rl
arl=c()
for(i in 1:l){
  arl=mean(rl[,i])
}
arl


















arl

我已经解决了，怎么弄回那1000个币呢

蓅☆煋語ぶ跡 发表于 2015-5-9 16:17
我已经解决了，怎么弄回那1000个币呢

重新编辑一下就可以了

oliyiyi 发表于 2015-5-9 22:47
重新编辑一下就可以了

就是问题再写一次？我不懂